{"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"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"***Next things to examine:***\n\n>\n> * Find out better ways to implement local cv scoring function\n> * Trying to implement logistic regression with neural networks (in order to try to tune it better)\n> * Explore other hyperparameter tuning techniques (GridSearchCV, RandomizedGridSearchCV, etc)\n> * Plot histograms for each model predictions. \n> * Feature Engineering. \n> * Ensembling: Blending, VotingClassifier, etc\n\n# <b>1 <span style='color:#3f4d63'>|</span> Introduction</b>\n\n<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>1.1 | Table of contents</b></p>\n</div>\n\n* **<span style = 'color:red'>Exploratory Data Analysis</span>**\n    * Basic information\n    * Statistical Analysis        \n    * Features' impact on failure    \n\n\n* **<span style = 'color:red'>Data Preprocessing</span>**    \n    * Encoding\n    * Scaling\n    * Feature Engineering\n        * Local CV Scoring Function\n        * PCA for FE\n    * Feature Selection (still working on)\n        * Mutual Information\n    \n    \n* **<span style = 'color:red'>Modeling</span>**    \n    * Using ML models\n        * PyCaret\n        * Hyperparameter Tuning. Optuna.\n            * Logistic Regression\n            * Linear Discriminant Analysis\n    * Neural Network. Keras approach\n    \n    \n* **<span style = 'color:red'>References</span>**\n\n<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>1.5 | Importing Packages</b></p>\n</div>","metadata":{}},{"cell_type":"code","source":"from IPython.display import clear_output\nimport os\nimport warnings\nfrom pathlib import Path\n\n# Basic libraries\nimport matplotlib.pyplot as plt\nimport numpy as np\nimport pandas as pd\n#import pandas_profiling as pp\nimport seaborn as sns\nimport scipy as sc\nfrom scipy import stats\n\n# Principal Component Analysis (PCA)\nfrom sklearn.decomposition import PCA\n\n#Mutual Information\nfrom sklearn.feature_selection import mutual_info_regression\n\n# Cross Validation\nfrom sklearn.model_selection import KFold, cross_val_score, StratifiedKFold, learning_curve, train_test_split\n\n# Plotly\nimport plotly.express as px\nfrom plotly.subplots import make_subplots\nimport plotly.figure_factory as ff\nimport plotly.offline as offline\nimport plotly.graph_objs as go\n\nwarnings.filterwarnings('ignore')\n\ndef load_data():\n    data_dir = Path(\"../input/tabular-playground-series-aug-2022\")\n    train = pd.read_csv(data_dir / \"train.csv\", index_col = [0])\n    test = pd.read_csv(data_dir / \"test.csv\", index_col = [0])\n    sample_submission = pd.read_csv(data_dir / 'sample_submission.csv')\n    return train, test, sample_submission\n\n\n# Function taken from Schubert de Abreu's work: \"Enthusiast to Data Professional - What Changes ?\"\n\"\"\"\nfunction annotation_helper(...)\n\nHelper for annotations in plotly. While reducing the amount of code to create an annotation, it also:\n- Allows us to provide the text into an array of \n  strings(one for each new line) instead of one really long <br> separated text param\n- Provides basic functionality for individual line spacing(s) between each line\n- Custom annotation rectangle\n- Basic debugging for annotation positioning\n\"\"\"\n\ndef annotation_helper(fig, texts, x, y, line_spacing, align=\"left\", bgcolor=\"rgba(0,0,0,0)\", borderpad=0, ref=\"axes\", xref=\"x\", yref=\"y\", width=100, debug = False):\n    \n    is_line_spacing_list = isinstance(line_spacing, list)\n    total_spacing = 0\n    \n    for index, text in enumerate(texts):\n        if is_line_spacing_list and index!= len(line_spacing):\n            current_line_spacing = line_spacing[index]\n        elif not is_line_spacing_list:\n            current_line_spacing = line_spacing\n        \n        fig.add_annotation(dict(\n            x= x,\n            y= y - total_spacing,\n            width = width,\n            showarrow=False,\n            text= text,\n            align= align,\n            borderpad=4 if debug == False else 0, # doesn't work with new background box implementation :S\n            xref= \"paper\" if ref==\"paper\" else xref,\n            yref= \"paper\" if ref==\"paper\" else yref,\n            \n            bordercolor= \"#222\",\n            borderwidth= 2 if debug == True else 0 # shows the actual borders of the annotation box\n        ))\n        \n        total_spacing  += current_line_spacing\n    \n    if bgcolor != \"rgba(0,0,0,0)\":\n        fig.add_shape(type=\"rect\",\n            xref= \"paper\" if ref==\"paper\" else xref,\n            yref= \"paper\" if ref==\"paper\" else yref,\n            xanchor = x, xsizemode = \"pixel\", \n            x0=-width/2, x1= +width/2, y0=y + line_spacing[-1], y1=y -total_spacing,\n            fillcolor= bgcolor,\n            line = dict(width=0))  \n      \n    if debug == True:\n        handle_annot_debug(fig, x, y, ref)\n        \ndef plot_feature_importance(importance,names,model_type):\n    \n    #Create arrays from feature importance and feature names\n    feature_importance = np.array(importance)\n    feature_names = np.array(names)\n    \n    #Create a DataFrame using a Dictionary\n    data={'feature_names':feature_names,'feature_importance':feature_importance}\n    fi_df = pd.DataFrame(data)\n    \n    #Sort the DataFrame in order decreasing feature importance\n    fi_df.sort_values(by=['feature_importance'], ascending=False,inplace=True)\n    fi_df = fi_df[fi_df.feature_importance > 0]\n    fig = px.bar(fi_df, x='feature_names', y='feature_importance', color=\"feature_importance\",\n             color_continuous_scale='Teal')\n    # General Styling\n    fig.update_layout(height=400, bargap=0.2,\n                  margin=dict(b=0),\n                  plot_bgcolor='rgb(242,242,242)',\n                  title = \"<span style='font-size:36px; font-family:Times New Roman'>Feature Importance Analysis</span>\",                  \n                  #paper_bgcolor = 'rgb(242,242,242)',\n                  font=dict(family=\"Times New Roman\", size= 14),\n                  hoverlabel=dict(font_color=\"floralwhite\"),\n                  showlegend=True)\n\n    fig.show()\n\ntrain, test, sample_submission = load_data()\nclear_output()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:01:02.967571Z","iopub.execute_input":"2022-08-07T00:01:02.968151Z","iopub.status.idle":"2022-08-07T00:01:03.201413Z","shell.execute_reply.started":"2022-08-07T00:01:02.968104Z","shell.execute_reply":"2022-08-07T00:01:03.200561Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# <b>2 <span style='color:#3f4d63'>|</span> Exploratory Data Analysis</b>\n\n<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>2.1 | Basic information</b></p>\n</div>\n\nTo start with, we're gonna take a brief view on the dataset given in order to get some basic information about it: \n\n* We're gonna show some sample rows of the dataframe.\n* Examine which type of features we've given\n* Find out some statistical information about each feature, as well as examining whether there are missing values.\n* Analyse how balanced is `failure`, our target feature. ","metadata":{}},{"cell_type":"code","source":"# Shape \nprint('Train set shape:', train.shape)\nprint('Test set shape:', test.shape)\nprint('\\n',train.dtypes, '\\n')\n\n# Let's first take a brief look onto the data we're given\ntrain.head()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:01:03.202917Z","iopub.execute_input":"2022-08-07T00:01:03.203314Z","iopub.status.idle":"2022-08-07T00:01:03.235517Z","shell.execute_reply.started":"2022-08-07T00:01:03.203272Z","shell.execute_reply":"2022-08-07T00:01:03.234547Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train.describe().style.background_gradient(cmap='Blues', axis = 1)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:03.237072Z","iopub.execute_input":"2022-08-07T00:01:03.237595Z","iopub.status.idle":"2022-08-07T00:01:03.369547Z","shell.execute_reply.started":"2022-08-07T00:01:03.237563Z","shell.execute_reply":"2022-08-07T00:01:03.368302Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test.describe().style.background_gradient(cmap='Blues', axis = 1)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:03.428696Z","iopub.execute_input":"2022-08-07T00:01:03.429109Z","iopub.status.idle":"2022-08-07T00:01:03.549560Z","shell.execute_reply.started":"2022-08-07T00:01:03.429079Z","shell.execute_reply":"2022-08-07T00:01:03.548251Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"***Insights:***\n\n* Taking a look at `Count` column above we appreciate that there are some missing values in both datasets. Let's examine it better.\n* Most of the numeric features tend to have a similar scale of values. However, `loading` and `measurement_17` break with this trendency. \n\n#### **Missing values**","metadata":{}},{"cell_type":"code","source":"# Reference for output: https://www.kaggle.com/code/ambrosm/tpsaug22-eda-which-makes-sense\ncols = test.columns\nprint('Missing values per feature: ')\npd.concat([train[cols].isna().sum().rename('Train'), test[cols].isna().sum().rename('Test')], axis=1).T.style.background_gradient(cmap='Blues', axis = 1)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:03.607622Z","iopub.execute_input":"2022-08-07T00:01:03.608026Z","iopub.status.idle":"2022-08-07T00:01:03.680539Z","shell.execute_reply.started":"2022-08-07T00:01:03.607995Z","shell.execute_reply":"2022-08-07T00:01:03.679131Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We should ask the following question: **Does these missing values mean something ?** In other words, as measurements features are related to various lab testing methods, a missing value could mean that the product would have failed already in a previous testing method. In fact, notice that as the later the measurement feature appears, the greater the amount of NaN values is: ","metadata":{}},{"cell_type":"code","source":"train[train.measurement_3.isnull() == True].loc[:, 'measurement_3':].head()","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:03.729925Z","iopub.execute_input":"2022-08-07T00:01:03.730570Z","iopub.status.idle":"2022-08-07T00:01:03.757391Z","shell.execute_reply.started":"2022-08-07T00:01:03.730537Z","shell.execute_reply":"2022-08-07T00:01:03.756394Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As we can appreciate, our conclusions does not fit properly with the data. We can observe a few samples where `measurement_3` has no data recorded, but most of the measurements already have it. Let's examine if having a missing value could have an impact on product's failure.","metadata":{}},{"cell_type":"code","source":"size = train[(train.measurement_3.isnull() == True)].shape[0]\nprint('Amount of samples with NaN value in measurement_3: ', size)\nprint('Failure %: ', train[(train.measurement_3.isnull() == True) & (train.failure == 1)].shape[0] / size)\nprint('\\nSamples with data recorded for measurement_3: ')\nprint('Failure %: ', train[(train.measurement_3.isnull() == False) & (train.failure == 1)].shape[0] / train.shape[0])","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:03.874498Z","iopub.execute_input":"2022-08-07T00:01:03.875124Z","iopub.status.idle":"2022-08-07T00:01:03.887527Z","shell.execute_reply.started":"2022-08-07T00:01:03.875090Z","shell.execute_reply":"2022-08-07T00:01:03.886586Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As we can notice, probability has change quite a bit. Therefore, let's examine this for the rest of the measurements features: ","metadata":{}},{"cell_type":"code","source":"nan_cols = [col for col in train.columns if train[col].isnull().sum() > 0]\nfor i, col in enumerate(nan_cols):\n    size = train[(train[col].isnull() == True)].shape[0]\n    print('Amount of samples with NaN value in {}: '.format(col), size)\n    print('Failure %: ', train[(train[col].isnull() == True) & (train.failure == 1)].shape[0] / size)\n    print('Samples with data recorded for {}: '.format(col))\n    print('Failure %: ', train[(train[col].isnull() == False) & (train.failure == 1)].shape[0] / train.shape[0])\n    print('\\n')","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:03.911434Z","iopub.execute_input":"2022-08-07T00:01:03.911889Z","iopub.status.idle":"2022-08-07T00:01:04.000152Z","shell.execute_reply.started":"2022-08-07T00:01:03.911838Z","shell.execute_reply":"2022-08-07T00:01:03.998639Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### **Filling Missing Values**","metadata":{}},{"cell_type":"code","source":"# Code taken from: https://www.kaggle.com/code/pourchot/keras-optuna-for-logistic-reg\n# Author: @pourchot\n# I'll examine this by my own when finishing the baseline DNN approach\n\nfrom sklearn.impute import KNNImputer\n\nna_col = [col for col in train.columns if train[col].isnull().sum() !=0]\nimputer = KNNImputer(n_neighbors = 3)\n\nfor df in [train, test]:\n    df[na_col] = imputer.fit_transform(df[na_col])\n    df.isnull().sum().sum()","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:04.004602Z","iopub.execute_input":"2022-08-07T00:01:04.005530Z","iopub.status.idle":"2022-08-07T00:01:57.481316Z","shell.execute_reply.started":"2022-08-07T00:01:04.005479Z","shell.execute_reply":"2022-08-07T00:01:57.479220Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### **Target's Balance**\n\nAs we can apprecite below our target feature is quite unbalanced. ","metadata":{}},{"cell_type":"code","source":"sns.countplot(train['failure'])","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:57.484624Z","iopub.execute_input":"2022-08-07T00:01:57.486613Z","iopub.status.idle":"2022-08-07T00:01:57.736564Z","shell.execute_reply.started":"2022-08-07T00:01:57.486553Z","shell.execute_reply":"2022-08-07T00:01:57.735055Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>2.2 | Statistical Analysis</b></p>\n</div>\n\n---\n### ***Correlations. Heatmap Chart***\n\n---\n\nIn this sub-section, our main aim will be to analyse the different relationships between each of the features. Due to it, we'll start by calculating their correlation coefficients and showing them in a heatmap chart. Thus, we'll be able to determine which features are linearly related.","metadata":{}},{"cell_type":"code","source":"corr = train.select_dtypes(['int','float']).corr()\ns = corr.unstack()\nso = s[s < 1.0].sort_values(ascending = False)\nprint('Most correlated features: ')\nmost_corr_f = ['measurement_17','measurement_8','measurement_5','attribute_3','attribute_2','measurement_1']\nso[abs(so) >= 0.25]","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:01:57.738198Z","iopub.execute_input":"2022-08-07T00:01:57.738559Z","iopub.status.idle":"2022-08-07T00:01:57.806066Z","shell.execute_reply.started":"2022-08-07T00:01:57.738527Z","shell.execute_reply":"2022-08-07T00:01:57.804733Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Let's now examine the whole correlation matrix by plotting in on a heatmap chart. ","metadata":{}},{"cell_type":"code","source":"fig, axes = plt.subplots(nrows = 2, ncols = 2, figsize = (22,16))\nfor i, df in enumerate([train, test]):\n    corr= df.select_dtypes(['int','float']).corr()\n    # Getting the Upper Triangle of the co-relation matrix\n    matrix = np.triu(corr)\n    \n    # Title\n    title = 'Train' if i == 0 else 'Test'\n    \n    # Heatmap without absolute values\n    sns.heatmap(corr, mask=matrix, center = 0, cmap = 'vlag', ax = axes[i][0]).set_title('{}Set without absolute values'.format(title))\n    # Heatmap with absolute values\n    sns.heatmap(abs(corr), mask=matrix, center = 0, cmap = 'vlag', ax = axes[i][1]).set_title('{}Set with absolute values'.format(title))\n\nfig.tight_layout(h_pad=1.2, w_pad=0.5)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:01:57.808950Z","iopub.execute_input":"2022-08-07T00:01:57.810245Z","iopub.status.idle":"2022-08-07T00:02:00.836651Z","shell.execute_reply.started":"2022-08-07T00:01:57.810201Z","shell.execute_reply":"2022-08-07T00:02:00.835357Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---\n### ***Continuous Features***\n\n---\n\nThe procedure will be the following: \n\n* Examine **how** continuous **data is distributed**. \n    * Distribution Charts. \n    * Statistical Non-Parametric Inference Techniques. ","metadata":{}},{"cell_type":"code","source":"from scipy import stats\n\nfigure = plt.figure(figsize = (22,12))\nfor j, df in enumerate([train.drop('failure',axis=1), test]):\n    label_name = 'train' if j == 0 else 'test'    \n    for i, col in enumerate(df.select_dtypes('float').columns): \n        plt.subplot(4, 4, i+1)\n        sns.distplot(df[col], fit=stats.norm, label = label_name)        \n        plt.legend()\n    \nfigure.tight_layout(h_pad=1.0, w_pad=0.5)\nplt.suptitle('Train and Test Distribution Plots', y=1.02)    \nplt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:00.838622Z","iopub.execute_input":"2022-08-07T00:02:00.839363Z","iopub.status.idle":"2022-08-07T00:02:15.504416Z","shell.execute_reply.started":"2022-08-07T00:02:00.839319Z","shell.execute_reply":"2022-08-07T00:02:15.503092Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Some of the features seem to be normally distributed. Let's check it. To do this, we're gonna use some **statistical non-parametric inference technique called Shapiro-Wilk Test**. This test is used to test whether a dataset is distributed normally or not. \n\n***Shapiro-Wilk Test***\n\nThe **null hypothesis** is that a sample $x_1\\hspace{0.1cm},\\hspace{0.1cm}\\cdots\\hspace{0.1cm},\\hspace{0.1cm}x_n$ comes from a normally distributed population. It was published in 1965 by Samuel Shapiro and Martin Wilk. **It is considered one of the most powerful tests for normality testing.** The test stadistic will be: \n\n$$W = \\frac{(\\sum_{i=1}^{n}a_{i}x_i)^2}{\\sum_{i=1}^{n}(x_i - \\bar{x})^2}$$\n\nwhere\n\n* $x_i$ is the number occupying the i-th position in the sample (with the sample ordered from smallest to largest).\n* $\\bar{x}$ is the sample mean. \n* Variables $a_i$ are calculated this way: \n\n$$(a_1, ... , a_n) = \\frac{m^T V^{-1}}{(m^T V^{-1}V^{-1}m)^{1/2}} \\hspace{2cm}m = (m_1 , ... , m_n)$$\n\nwhere $m_1 , ... , m_n$ are the mean values of the ordered statistic, of independent and identically distributed random variables, sampled from normal distributions and $V$ denotes the covariance matrix of that order statistic. **The null hypothesis is rejected if W is too small. The value of W can range from 0 to 1.**","metadata":{}},{"cell_type":"code","source":"from scipy.stats import shapiro\nfrom termcolor import colored\n\n# References: https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.shapiro.html\n#             https://www.kaggle.com/code/kartushovdanil/tps-jul-22-advanced-2-sol (Author: Torch me)\n\ndef shapiro_wilk_test(type = 'float'): \n    non_normal_f = []\n    for i, col in enumerate(train.select_dtypes(type).columns): \n        test_statistic, p_value = shapiro(train[col])  \n        alpha = 0.05\n        if p_value > alpha: \n            result = colored('Accepted', 'green')  \n        else:\n            result = colored('Rejected','red') \n            non_normal_f.append(col)\n        print('Feature: {}\\t Hypothesis: {}'.format(col, result))\n    return non_normal_f\n\nnon_normal_f = shapiro_wilk_test()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:15.506315Z","iopub.execute_input":"2022-08-07T00:02:15.507515Z","iopub.status.idle":"2022-08-07T00:02:15.567924Z","shell.execute_reply.started":"2022-08-07T00:02:15.507467Z","shell.execute_reply":"2022-08-07T00:02:15.566632Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Only features from `measurement_3` to `measurement_9` are normally distributed. \n\n***Poisson Dispersion Test*** (Credits to: [TPS July22. Author: @samuelcortinhas](https://www.kaggle.com/code/javigallego/tps-july-22-unsupervised-clustering-949a64/edit))\n\nThis test is used to determine whether a feature is distributed according to a **Poisson** distribution or not. The null hypothesis is that \n\n$$\nX_i \\sim Po(\\lambda) \\quad \\text{for every i=1, $\\ldots$, n}\n$$\n\nThis is the **most common** test used for verifying a Poisson distribution. The test statistic (called dispersion) is\n\n$$\nD = \\sum_{i=1}^{n} \\frac{(X_i - \\bar{X})^2}{\\bar{X}},\n$$\n\nwhere\n\n* $X_i$ is the number from the i-th sample point (order doesn't matter)\n* $\\bar{X}$ is the sample mean.\n\nNote that the **expected value** of this statistic is $\\mathbb{E}(D) = \\frac{(n-1) Var(X_i)}{E(X_i)} = \\frac{(n-1) \\lambda}{\\lambda}  = n-1$, since the mean and variance of a Poisson distribution is the rate $\\lambda$. If $D$ is too '**far away**' from the expected value of $n-1$, then we **reject** the null hypothesis. \n\nMore formally, $D$ has a **chi-squared** distribution with $n-1$ **degrees of freedom** under the null hypothesis. We determine the **critical values** by using a **two-tailed** test with significance level $\\alpha=5\\%$.","metadata":{}},{"cell_type":"code","source":"# Credits: https://www.kaggle.com/code/javigallego/tps-july-22-unsupervised-clustering-949a64/edit\n#          Author: @samuelcortinhas\n\nfrom scipy.stats import chi2\nfrom scipy.stats import poisson\n\ndef poisson_dispersion_test():\n    for col in non_normal_f:\n        # Parameters\n        alpha = 0.05                  # significance level\n        n = len(train[col])            # sample size\n        df = n-1                      # degrees of freedom\n\n        # Statistics\n        mu = train[col].mean()               # sample mean\n        D = ((train[col]-mu)**2).sum()/mu    # test statistic\n\n        # Two-tailed test\n        q_lower = alpha/2\n        q_upper = (1-alpha)/2\n\n        # percentile point function = inverse of cdf\n        chi2_crit_lower = chi2.ppf(q_lower, df)\n        chi2_crit_upper = chi2.ppf(q_upper, df)\n\n        if (D<chi2_crit_lower) or (D>chi2_crit_upper):\n            result = colored('Rejected', 'red')\n        else:\n            result = colored('Accepted', 'green')\n        print('Feature: {}\\t Hypothesis: {}'.format(col, result))\n        #print('D:',int(D),', chi2_crit_lower:',int(chi2_crit_lower),', chi2_crit_upper:',int(chi2_crit_upper),'\\n')\n    \npoisson_dispersion_test()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:15.569911Z","iopub.execute_input":"2022-08-07T00:02:15.570712Z","iopub.status.idle":"2022-08-07T00:02:15.597736Z","shell.execute_reply.started":"2022-08-07T00:02:15.570653Z","shell.execute_reply":"2022-08-07T00:02:15.595439Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---\n### ***Discrete Features***\n\n---\n\nThe preocedure will be the following. We'll start by analysing `product_code`. This feature distinguish one type of product from another. However, we'll notice that to dig into this feature's sense, we'll need to analyse `attribute` features. To finish with, we'll examine `measurement` integer features. \n\n#### **Examining each type of product**\n\nThe procedure to analyse this feature will be as follows: \n\n* Examining both training and testing set values for this feature. \n* Finding out its impact on `failure`. \n* Digging into this feature's sense. \n    * Analyse `attribute` features. ","metadata":{}},{"cell_type":"code","source":"pd_codes = []\nfor df in [train, test]:\n    pd_codes.append(df.groupby('product_code').count().loc[:, df.columns[1]] / df.shape[0])\n    \npd_codes = pd.DataFrame(pd_codes).transpose()\npd_codes.columns = ['train','test']\n\n\nfig, axes = plt.subplots(nrows = 1, ncols = 2, sharey = True, figsize = (16,4))\n\nfor i, col in enumerate(['train','test']):\n    sns.barplot(x = pd_codes.index, y = pd_codes[col], ax = axes[i])\n\nfigure.tight_layout(h_pad=1.0, w_pad=0.5)\nplt.suptitle('Product_Code Percentage', y=1.02)\nfig.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:15.601997Z","iopub.execute_input":"2022-08-07T00:02:15.603644Z","iopub.status.idle":"2022-08-07T00:02:16.969938Z","shell.execute_reply.started":"2022-08-07T00:02:15.603542Z","shell.execute_reply":"2022-08-07T00:02:16.968480Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Insights**\n* We do not have either of the `product_code` of the training set in the testing one. \n* Nine different types of products. \n* **We may not use this feature for modeling !**","metadata":{}},{"cell_type":"code","source":"failures = []\nfor i, val in enumerate(train.product_code.unique()): \n    for j in [0,1]: \n        mask = (train.product_code == val) & (train.failure == j)\n        failure_percentage = train[mask].shape[0] / train[train.product_code == val].shape[0]\n        failures.append([val, j,failure_percentage])\n    \nfailures = pd.DataFrame(failures, columns = ['product_code','Failure','Percentage'])\nsns.barplot(data = failures, x = 'product_code', y = 'Percentage', hue = 'Failure')","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:16.972148Z","iopub.execute_input":"2022-08-07T00:02:16.972914Z","iopub.status.idle":"2022-08-07T00:02:17.327609Z","shell.execute_reply.started":"2022-08-07T00:02:16.972869Z","shell.execute_reply":"2022-08-07T00:02:17.326735Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We appreciate a similar tendency as noticed early. Each `product_type` has a failure probability around 0.2. We were aiming to notice whether there are some types of products which are more keen to fail. However, it's not the case. \n\n##### **Correlations per product code**","metadata":{}},{"cell_type":"code","source":"fig, axes = plt.subplots(nrows = 3, ncols = 3, figsize = (22,18))\n\nfor i, val in enumerate(train.product_code.unique()): \n    corr = train[train.product_code == val].corr()\n    matrix = np.triu(corr)\n    sns.heatmap(corr, mask=matrix, center = 0, cmap = 'vlag', ax = axes[i // 3][i%3]).set_title('Train Product Code: {}'.format(val))\n    \nfor val in test.product_code.unique(): \n    i += 1\n    corr = test[test.product_code == val].corr()\n    matrix = np.triu(corr)\n    sns.heatmap(corr, mask=matrix, center = 0, cmap = 'vlag', ax = axes[i // 3][i%3]).set_title('Test Product Code: {}'.format(val))    \n\nfig.tight_layout(h_pad=1.2, w_pad=0.5)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:17.331701Z","iopub.execute_input":"2022-08-07T00:02:17.332280Z","iopub.status.idle":"2022-08-07T00:02:25.745393Z","shell.execute_reply.started":"2022-08-07T00:02:17.332247Z","shell.execute_reply":"2022-08-07T00:02:25.744118Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### **Digging into product_code Sense**\n\nLet's dig a bit more into this feature's sense. If we recap, we have two `object` features. Both `attribute_0` and `attribute_1`. I tend to believe that each product code represents a concrete combination of this two features. In fact, it could involve as well the rest of the discrete features. But we'll focus on this later: ","metadata":{}},{"cell_type":"code","source":"def attribute_info(col): \n    print('All possible values for {}:'.format(col))\n    for j, df in enumerate([train, test]):\n        name = 'Train' if j == 0 else 'Test '\n        print('{}: {}'.format(name, df[col].unique()))\n        \n    print('\\n')\n    fig, axes = plt.subplots(nrows = 1, ncols = 2, figsize = (16,4))\n    for i, df in enumerate([train, test]):     \n        sns.countplot(df[col], ax = axes[i])        \n        \nattribute_info('attribute_0')","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:25.746995Z","iopub.execute_input":"2022-08-07T00:02:25.747493Z","iopub.status.idle":"2022-08-07T00:02:26.073277Z","shell.execute_reply.started":"2022-08-07T00:02:25.747449Z","shell.execute_reply":"2022-08-07T00:02:26.071917Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"attribute_info('attribute_1')","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:26.074858Z","iopub.execute_input":"2022-08-07T00:02:26.075331Z","iopub.status.idle":"2022-08-07T00:02:26.387776Z","shell.execute_reply.started":"2022-08-07T00:02:26.075298Z","shell.execute_reply":"2022-08-07T00:02:26.386513Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As you may have noticed we only have 4 types of material, gathering together both features. From my point of view, this could mean that each of the products are made with two of these materials. This could make sense for creating the `product_code` feature. Let's examine this:  ","metadata":{}},{"cell_type":"code","source":"for df in [train, test]:\n    for val in df.product_code.unique():\n        print('Product code: {}'.format(val))\n        print(df[df.product_code == val]['attribute_0'].value_counts())\n        print(df[df.product_code == val]['attribute_1'].value_counts())\n        print('\\t')","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:26.389198Z","iopub.execute_input":"2022-08-07T00:02:26.389606Z","iopub.status.idle":"2022-08-07T00:02:26.472924Z","shell.execute_reply.started":"2022-08-07T00:02:26.389576Z","shell.execute_reply":"2022-08-07T00:02:26.471552Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We may conclude that our intuitions were right. However, you may have noticed one interesting aspect. Denoting `attribute_0` and `attribute_1` with **layer1** and **layer2** respectively, **<span style='color:red'>there are some types sharing the same mixture for its layers</span>**. \n\n* **What does this mean ?** There has to be another criterion to distinguish this types. \n* **Which criterions ?** The ***rest of the attribute features***.\n\nLet's take a brief view on them: ","metadata":{}},{"cell_type":"code","source":"attribute_info('attribute_2')","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:26.474692Z","iopub.execute_input":"2022-08-07T00:02:26.475189Z","iopub.status.idle":"2022-08-07T00:02:26.777747Z","shell.execute_reply.started":"2022-08-07T00:02:26.475147Z","shell.execute_reply":"2022-08-07T00:02:26.776716Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"attribute_info('attribute_3')","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:26.779270Z","iopub.execute_input":"2022-08-07T00:02:26.779608Z","iopub.status.idle":"2022-08-07T00:02:27.071811Z","shell.execute_reply.started":"2022-08-07T00:02:26.779579Z","shell.execute_reply":"2022-08-07T00:02:27.070597Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Thus, proceeding in the same manner as we did before we'll find out how to distinguish one type from another: \n\n> ***Reminder:*** *pairs D-I and F-G are the ones having the same mixture of materials for `attribute_0` and `attribute_1`. In other words, we have not been able to make a distinction between them yet.*","metadata":{}},{"cell_type":"code","source":"for df in [train, test]:\n    for val in df.product_code.unique():\n        print('Product code: {}'.format(val))\n        print(df[df.product_code == val]['attribute_2'].value_counts())\n        print(df[df.product_code == val]['attribute_3'].value_counts())\n        print('\\t')","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:27.073214Z","iopub.execute_input":"2022-08-07T00:02:27.073554Z","iopub.status.idle":"2022-08-07T00:02:27.153523Z","shell.execute_reply.started":"2022-08-07T00:02:27.073523Z","shell.execute_reply":"2022-08-07T00:02:27.152134Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now we are able to distinguish properly each of the types of products from the others. In fact, for those having the same mixture of materials (products of type D-I and F-G) `attribute_2` and `attribute_3` make finally the distinction. \n\n#### **Measurement integer features**\n\nTo finish with let's make a brief analysis on this last kind of features: ","metadata":{}},{"cell_type":"code","source":"def integer_measurements_histogram(type, r = 4, height = 10): \n    figure, axes = plt.subplots(nrows = r, ncols = 2, figsize = (22,height))\n    for j, df in enumerate([train.drop(['product_code','failure'],axis=1), test.drop('product_code',axis=1)]):\n        label_name = 'Train' if j == 0 else 'Test'    \n        for i, col in enumerate(['measurement_{}'.format(i) for i in range(3)]): \n            sns.countplot(df[col], ax = axes[i][j])        \n            axes[0][j].set_title('{}'.format(label_name))                \n\n    figure.tight_layout(h_pad=1.0, w_pad=0.5)\n    plt.suptitle('{} features Plots'.format(type), y=1.02)    \n    plt.show()\n    \ninteger_measurements_histogram('int', 3)    ","metadata":{"_kg_hide-input":false,"execution":{"iopub.status.busy":"2022-08-07T00:02:27.156233Z","iopub.execute_input":"2022-08-07T00:02:27.157167Z","iopub.status.idle":"2022-08-07T00:02:29.410823Z","shell.execute_reply.started":"2022-08-07T00:02:27.157117Z","shell.execute_reply":"2022-08-07T00:02:29.409468Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"To finish this sub-section let's find out if some features' distribution is a Poisson Distribution (discrete measurements seem to have it). ","metadata":{}},{"cell_type":"code","source":"print('Shapiro-Wilk Test: ')\nnon_normal_f = shapiro_wilk_test('int')\nprint('\\nPoisson Dispersion Test: ')\npoisson_dispersion_test()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:29.412670Z","iopub.execute_input":"2022-08-07T00:02:29.413888Z","iopub.status.idle":"2022-08-07T00:02:29.439025Z","shell.execute_reply.started":"2022-08-07T00:02:29.413844Z","shell.execute_reply":"2022-08-07T00:02:29.437816Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---\n### **Skewness and Kurtosis**\n---\n\n> [Skewness and Kurtosis Tutorial](https://www.universoformulas.com/estadistica/descriptiva/asimetria-curtosis/)\n\n**Skewness**","metadata":{}},{"cell_type":"code","source":"from scipy.stats import skew, kurtosis\n\nfig = make_subplots(rows = 1, cols = 1)\nfor i, df in enumerate([train.drop('failure',axis = 1), test]):\n\n    color_scale = px.colors.qualitative.G10 if i == 0 else px.colors.qualitative.T10\n    label_name = 'Train' if i == 0 else 'Test'\n    \n    quantitative = df.select_dtypes(['int','float']).columns\n    skew_features = pd.DataFrame(df[quantitative].apply(lambda x : skew(x))).reset_index()\n    skew_features.columns = ['features','skewness']\n    \n    fig.add_trace(go.Bar(x = skew_features['features'] ,y=skew_features['skewness'], \n                         marker = dict(color = color_scale[-1]), name = label_name), \n                  row = 1, col = 1)\n\n# General Styling\nfig.update_layout(height=400, bargap=0.2,\n                  margin=dict(b=0, t=120),\n                  plot_bgcolor='rgb(242,242,242)',\n                  #paper_bgcolor = 'rgb(242,242,242)',\n                  font=dict(family=\"Times New Roman\", size= 14),\n                  hoverlabel=dict(font_color=\"floralwhite\"),\n                  showlegend=True)\n\n#Title\ntext = [\n    \"<span style='font-size:35px; font-family:Times New Roman'>Skewness Chart</span>\"\n]\nannotation_helper(fig, text, -0.043, 1.5, [0.095,0.095,0.065],ref=\"paper\", width=1300)\n\n#Subtitle\ntext = [\n    \"<span style='font-size:14px; font-family:Helvetica'> Skewness is a measure of the asymmetry of the probability distribution of a real-valued random variable about its mean. The skewness value can be positive, zero, negative, or undefined.</span>\",\n    \"<span style='font-size:14px; font-family:Helvetica'> For a unimodal distribution, negative skew commonly indicates that the tail is on the left side of the distribution, and positive skew indicates that the tail is on the right</span>\",\n]\nannotation_helper(fig, text, -0.043, 1.3, [0.075,0.05,0.065],ref=\"paper\", width=1300)\n\nfig.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:29.441072Z","iopub.execute_input":"2022-08-07T00:02:29.441908Z","iopub.status.idle":"2022-08-07T00:02:29.761813Z","shell.execute_reply.started":"2022-08-07T00:02:29.441861Z","shell.execute_reply":"2022-08-07T00:02:29.760821Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Kurtosis**","metadata":{}},{"cell_type":"code","source":"from scipy.stats import skew, kurtosis\n\nfig = make_subplots(rows = 1, cols = 1)\nfor i, df in enumerate([train.drop('failure',axis = 1), test]):\n\n    color_scale = px.colors.qualitative.G10 if i == 0 else px.colors.qualitative.T10\n    label_name = 'Train' if i == 0 else 'Test'\n    \n    quantitative = df.select_dtypes(['int','float']).columns\n    skew_features = pd.DataFrame(df[quantitative].apply(lambda x : kurtosis(x))).reset_index()\n    skew_features.columns = ['features','kurtosis']\n    \n    fig.add_trace(go.Bar(x = skew_features['features'] ,y=skew_features['kurtosis'], \n                         marker = dict(color = color_scale[-1]), name = label_name), \n                  row = 1, col = 1)\n\n# General Styling\nfig.update_layout(height=400, bargap=0.2,\n                  margin=dict(b=0, t=120),\n                  plot_bgcolor='rgb(242,242,242)',\n                  #paper_bgcolor = 'rgb(242,242,242)',\n                  font=dict(family=\"Times New Roman\", size= 14),\n                  hoverlabel=dict(font_color=\"floralwhite\"),\n                  showlegend=True)\n\n#Title\ntext = [\n    \"<span style='font-size:35px; font-family:Times New Roman'>Kurtosis Chart</span>\"\n]\nannotation_helper(fig, text, -0.043, 1.5, [0.095,0.095,0.065],ref=\"paper\", width=1300)\n\n#Subtitle\ntext = [\n    \"<span style='font-size:14px; font-family:Helvetica'> Kurtosis is a measure of the tailedness of the probability distribution of a real-valued random variable. Like skewness, it describes the shape of a probability distribution and there are different</span>\",\n    \"<span style='font-size:14px; font-family:Helvetica'> ways of quantifying it for a theoretical distribution and corresponding ways of estimating it from a sample from a population. Different measures of kurtosis may have different interpretations.</span>\",\n]\nannotation_helper(fig, text, -0.043, 1.3, [0.075,0.05,0.065],ref=\"paper\", width=1300)\n\nfig.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:29.763269Z","iopub.execute_input":"2022-08-07T00:02:29.763806Z","iopub.status.idle":"2022-08-07T00:02:29.828736Z","shell.execute_reply.started":"2022-08-07T00:02:29.763774Z","shell.execute_reply":"2022-08-07T00:02:29.827536Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>2.3 | Values' impact on Failure </b></p>\n</div>\n\nLet's now examine whether features' values have an impact on failure probability. To do so, we're going to plot failure mean per each of the groups for a single feature. But first, let's examine distribution of each numerical feature making a distinction between failure values. ","metadata":{}},{"cell_type":"code","source":"figure = plt.figure(figsize = (22,12))\nfor i, col in enumerate(train.drop('failure', axis = 1).select_dtypes(['int','float']).columns): \n    plt.subplot(6, 4, i+1)\n    sns.kdeplot(train[col], hue = train['failure'])        \n    \nfigure.tight_layout(h_pad=1.0, w_pad=0.5)\nplt.suptitle('Distribution Distinguishing Failures', y=1.02)    \nplt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:29.830615Z","iopub.execute_input":"2022-08-07T00:02:29.831003Z","iopub.status.idle":"2022-08-07T00:02:36.774648Z","shell.execute_reply.started":"2022-08-07T00:02:29.830974Z","shell.execute_reply":"2022-08-07T00:02:36.773081Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"From this previous chart, we can appreciate that initially numerical features do not really help us a lot distinguishing whether a product is gonna suffer a failure or not. ","metadata":{}},{"cell_type":"code","source":"figure, axes = plt.subplots(nrows = 2, ncols = 4, figsize = (22,10), sharey = True)\n\nfor i, col in enumerate(train.drop('failure',axis =1).select_dtypes(['int','object']).columns): \n    group = train[[col,'failure']].groupby(col).mean()    \n    sns.lineplot(x = group.index, y = group['failure'], ax = axes[i // 4][i % 4])        \n    plt.ylim(0,1)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:36.776587Z","iopub.execute_input":"2022-08-07T00:02:36.776946Z","iopub.status.idle":"2022-08-07T00:02:37.834582Z","shell.execute_reply.started":"2022-08-07T00:02:36.776918Z","shell.execute_reply":"2022-08-07T00:02:37.833209Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As observed, for the first five features the probability of failure is almost the same. However, this changes when talking about the last ones. Interesting insights to notice: \n\n* Samples with `measurement_0` $\\geq$ 26 have 0.0 probability of failure. \n* Samples with `measurement_1` is equal to 26 or 28 have 0.0 probability of failure. On the other hand, if this value is 27 or 29 this probability increases to 1.0. \n* Samples with `measurement_2` = 22 have 1.0 probability of failure. On the other hand, if this value is 23 this probability decreases to 0.0. \n\nHowever, we must take into consideration that all of these values are almost at the extremes of the distribution. As distributions have a gaussian bell shape, it means that are non-common values for a sample. ","metadata":{}},{"cell_type":"code","source":"print('Extreme values percentage: ',train[train.measurement_0 >= 25].shape[0])\nprint('Extreme values percentage: ',train[train.measurement_1 >= 26].shape[0])\nprint('Extreme values percentage: ',train[train.measurement_2 >= 22].shape[0])","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:37.836473Z","iopub.execute_input":"2022-08-07T00:02:37.836939Z","iopub.status.idle":"2022-08-07T00:02:37.848992Z","shell.execute_reply.started":"2022-08-07T00:02:37.836895Z","shell.execute_reply":"2022-08-07T00:02:37.847609Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"In conclusion, if we do not take into account these unusual values we have the following: **failure probability oscilate between 0.2 for every feature and every possible value.** ","metadata":{}},{"cell_type":"code","source":"figure, axes = plt.subplots(nrows = 2, ncols = 4, figsize = (22,8), sharey = True)\n\nfor i, col in enumerate(train.drop(['product_code','failure'],axis =1).select_dtypes(['int','object']).columns): \n    if df[col].dtype == 'int': \n        mask = train[col] <= train[col].quantile(0.95)\n        group = train[mask][[col,'failure']].groupby(col).mean()    \n    else: \n        group = train[[col,'failure']].groupby(col).mean()\n    sns.lineplot(x = group.index, y = group['failure'], ax = axes[i // 4][i % 4])        \n    plt.ylim(0,1)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:37.850731Z","iopub.execute_input":"2022-08-07T00:02:37.851678Z","iopub.status.idle":"2022-08-07T00:02:38.889380Z","shell.execute_reply.started":"2022-08-07T00:02:37.851636Z","shell.execute_reply":"2022-08-07T00:02:38.888334Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---\n\n#### ***To be continued***\n\nFrom this section onwards I'll let just basic titles and explanations, apart from the code obviously.  I'll keep improving project's markdown and styling on next versions. \n\n---\n\n# <b>3 <span style='color:#3f4d63'>|</span> Data Preprocessing</b>\n\n<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>3.1 | Encoding </b></p>\n</div>","metadata":{}},{"cell_type":"code","source":"from sklearn.preprocessing import LabelEncoder\n\nfor col in test.select_dtypes('object').columns:\n    for df in [train, test]:\n        le = LabelEncoder()\n        df[col] = le.fit_transform(df[col])","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:38.890908Z","iopub.execute_input":"2022-08-07T00:02:38.891995Z","iopub.status.idle":"2022-08-07T00:02:38.939153Z","shell.execute_reply.started":"2022-08-07T00:02:38.891958Z","shell.execute_reply":"2022-08-07T00:02:38.937768Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>3.2 | Scaling </b></p>\n</div>","metadata":{}},{"cell_type":"code","source":"from sklearn.preprocessing import StandardScaler\n\nscaler = StandardScaler()\n\ntrain_scaled = pd.DataFrame(scaler.fit_transform(train.drop('failure',axis = 1))    )\ntrain_scaled.columns = train.columns[:-1]\ntrain_scaled = pd.concat([train_scaled, train['failure']], axis = 1)                                     \n\ntest_scaled = pd.DataFrame(scaler.fit_transform(test))\ntest_scaled.columns = test.columns ","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:38.940687Z","iopub.execute_input":"2022-08-07T00:02:38.941835Z","iopub.status.idle":"2022-08-07T00:02:38.986694Z","shell.execute_reply.started":"2022-08-07T00:02:38.941798Z","shell.execute_reply":"2022-08-07T00:02:38.985270Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>3.3 |  Feature Engineering</b></p>\n</div>\n\n> This is an initial version of this section. It'll be uploaded properly in next versions. \n\n---\n#### **Local CV Scoring Function**\n---\n\nFirst thing we're gonna implement is our local cv scoring function. We'll use it for examining how much impact new features created have on modeling. ","metadata":{}},{"cell_type":"code","source":"from sklearn.discriminant_analysis import LinearDiscriminantAnalysis as lda\nfrom sklearn.linear_model import LogisticRegression as lr\nfrom sklearn.ensemble import GradientBoostingClassifier as gbc\nfrom sklearn.metrics import roc_auc_score\nfrom sklearn.model_selection import train_test_split, StratifiedGroupKFold, GroupKFold\n\ndef lcv_scoring(X, y, model=lr()):\n    train_scores = []\n    val_scores = []        \n\n    kf = GroupKFold(n_splits=5)   # group by product_code\n    for fold, (idx_train, idx_valid) in enumerate(kf.split(X, y, groups = X['product_code'])):\n\n        X_train = X.drop('product_code',axis=1).iloc[idx_train]\n        X_valid = X.drop('product_code',axis=1).iloc[idx_valid]\n        y_train = y[idx_train]\n        y_valid = y[idx_valid]\n\n        model.fit(X = X_train, y = y_train)        \n        \n        y_train_pred = model.predict(X_train)\n        y_valid_pred = model.predict(X_valid)\n        train_auc_score = roc_auc_score(y_train, y_train_pred)\n        valid_auc_score = roc_auc_score(y_valid, y_valid_pred)                    \n        train_scores.append(train_auc_score)\n        val_scores.append(valid_auc_score)\n        \n    print(f'Average Accuracy Score of Train: {np.mean(train_scores)}')\n    print(f'Average Accuracy Score of Validation: {np.mean(val_scores)}')\n    return np.mean(val_scores) \n\n\ndef score_dataset():\n    X = train_scaled.drop('failure', axis = 1)\n    y = train_scaled.failure\n\n    baseline_score = lcv_scoring(X, y)\n    print(f\"ROC AUC Score: {baseline_score:.5f}\")\n    \nscore_dataset()","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:38.988911Z","iopub.execute_input":"2022-08-07T00:02:38.989415Z","iopub.status.idle":"2022-08-07T00:02:40.002652Z","shell.execute_reply.started":"2022-08-07T00:02:38.989370Z","shell.execute_reply":"2022-08-07T00:02:40.000751Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"---\n#### **PCA for Feature Engineering**\n---\n\nWe're going to study correlational structure of our data. To do so, we're going to use PCA. Let's consider only those features having the largest correlation coefficients.","metadata":{}},{"cell_type":"code","source":"from sklearn.preprocessing import *\nfrom sklearn.decomposition import PCA\n\ndef apply_pca(X, transformer = False, components = -1):\n    aux = X.copy()\n    if transformer:\n        X = pd.DataFrame(transformer.fit_transform(X))\n        X.columns = aux.columns    \n    # Create principal components\n    if components == -1:\n        pca = PCA()\n    else:\n        pca = PCA(n_components = components)\n        \n    X_pca = pca.fit_transform(X)\n    # Convert to dataframe\n    component_names = [f\"PC{i+1}\" for i in range(X_pca.shape[1])]\n    X_pca = pd.DataFrame(X_pca, columns=component_names)\n    # Create loadings\n    loadings = pd.DataFrame(\n        pca.components_.T,  # transpose the matrix of loadings\n        columns=component_names,  # so the columns are the principal components\n        index=X.columns,  # and the rows are the original features\n    )\n    return pca, X_pca, loadings\n\n\ndef plot_variance(pca, width=8, dpi=100):\n    # Create figure\n    fig, axs = plt.subplots(1, 2, sharey = True)\n    n = pca.n_components_\n    grid = np.arange(1, n + 1)\n    # Explained variance\n    evr = pca.explained_variance_ratio_\n    axs[0].bar(grid, evr)\n    axs[0].set(\n        xlabel=\"Component\", title=\"% Explained Variance\")    # ylim = (0.0, 1.0) o sino sharey en plt.subplots\n    cv = np.cumsum(evr)\n    axs[1].plot(np.r_[0, grid], np.r_[0, cv], \"o-\")\n    axs[1].set(\n        xlabel=\"Component\", title=\"% Cumulative Variance\"\n    )\n    # Set up figure\n    fig.set(figwidth=8, dpi=100)\n    return axs\n\nfeatures = []\nfeatures.append(\"attribute_2\")\nfeatures.append(\"attribute_3\")\nfeatures.append(\"measurement_8\")\nfeatures.append(\"measurement_17\")\nfeatures.append(\"measurement_6\")\n\nX = train_scaled.copy()\nX = X.loc[:, features]\n\npca, X_pca, loadings = apply_pca(X)\nplot_variance(pca)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:40.018300Z","iopub.execute_input":"2022-08-07T00:02:40.020328Z","iopub.status.idle":"2022-08-07T00:02:40.490317Z","shell.execute_reply.started":"2022-08-07T00:02:40.020252Z","shell.execute_reply":"2022-08-07T00:02:40.489078Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"loadings","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:40.492382Z","iopub.execute_input":"2022-08-07T00:02:40.492733Z","iopub.status.idle":"2022-08-07T00:02:40.507432Z","shell.execute_reply.started":"2022-08-07T00:02:40.492703Z","shell.execute_reply":"2022-08-07T00:02:40.506095Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_scaled['feature1'] = train_scaled['measurement_17'] * train_scaled['measurement_8']\ntrain_scaled['feature2'] = train_scaled['measurement_17'] + train_scaled['measurement_6']\nscore_dataset()","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:40.509631Z","iopub.execute_input":"2022-08-07T00:02:40.510293Z","iopub.status.idle":"2022-08-07T00:02:41.330724Z","shell.execute_reply.started":"2022-08-07T00:02:40.510126Z","shell.execute_reply":"2022-08-07T00:02:41.329397Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_scaled['feature3'] = train_scaled['measurement_10'] / train_scaled['measurement_12']\nscore_dataset()","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:41.332795Z","iopub.execute_input":"2022-08-07T00:02:41.333617Z","iopub.status.idle":"2022-08-07T00:02:43.707630Z","shell.execute_reply.started":"2022-08-07T00:02:41.333561Z","shell.execute_reply":"2022-08-07T00:02:43.706420Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_scaled['feature4'] = train_scaled['measurement_17'] * train_scaled['measurement_9']\nscore_dataset()","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:43.709520Z","iopub.execute_input":"2022-08-07T00:02:43.710225Z","iopub.status.idle":"2022-08-07T00:02:45.949503Z","shell.execute_reply.started":"2022-08-07T00:02:43.710181Z","shell.execute_reply":"2022-08-07T00:02:45.948274Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>3.4 |  Feature Selection</b></p>\n</div>\n\n> This is an initial version of this section. It'll be uploaded properly in next versions. \n\n---\n#### **Mutual Information**\n---","metadata":{}},{"cell_type":"code","source":"from sklearn.feature_selection import mutual_info_regression\n\ndef make_mi_scores(X, y):\n    X = X.copy()\n    for colname in X.select_dtypes([\"object\"]):\n        X[colname], _ = X[colname].factorize()\n    # All discrete features should now have integer dtypes\n    #discrete_features = [pd.api.types.is_integer_dtype(t) for t in X.dtypes]\n    mi_scores = mutual_info_regression(X, y, random_state=0)\n    mi_scores = pd.Series(mi_scores, name=\"MI Scores\", index=X.columns)\n    mi_scores = mi_scores.sort_values(ascending=False)\n    return mi_scores\n\ny = train_scaled['failure']\nx = train_scaled.drop('failure', axis=1)\nmi_scores = make_mi_scores(x, y)\nmi_scores = pd.DataFrame(mi_scores).reset_index().rename(columns={'index':'Feature'})\n\nfig = px.bar(mi_scores, x='MI Scores', y='Feature', color=\"MI Scores\",\n             color_continuous_scale='Blues')\n\n# General Styling\nfig.update_layout(height=750, bargap=0.2, xaxis={'categoryorder':'category ascending'},\n                  margin=dict(b=0, t=120),\n                  plot_bgcolor='rgb(242,242,242)',\n                  #paper_bgcolor = 'rgb(242,242,242)',\n                  font=dict(family=\"Times New Roman\", size= 14),\n                  hoverlabel=dict(font_color=\"floralwhite\"),\n                  showlegend=False)\n\n#Title\ntext = [\n    \"<span style='font-size:35px; font-family:Times New Roman'>Mutual Information</span>\"\n]\nannotation_helper(fig, text, -0.043, 1.225, [0.095,0.095,0.065],ref=\"paper\", width=1300)\n\n#Subtitle\ntext = [\n    \"<span style='font-size:14px; font-family:Helvetica'> Mutual information describes relationships in terms of uncertainty. The mutual information (MI) between two quantities is a measure of the extent to which knowledge of one quantity reduces uncertainty about the other. </span>\",\n    \"<span style='font-size:14px; font-family:Helvetica'> If you knew the value of a feature, how much more confident would you be about the target? </span>\",\n]\nannotation_helper(fig, text, -0.043, 1.135, [0.035,0.05,0.065],ref=\"paper\", width=1300)\n\nfig.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:02:45.951596Z","iopub.execute_input":"2022-08-07T00:02:45.952420Z","iopub.status.idle":"2022-08-07T00:02:51.847381Z","shell.execute_reply.started":"2022-08-07T00:02:45.952372Z","shell.execute_reply":"2022-08-07T00:02:51.845738Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_scaled['feature1'] = test_scaled['measurement_17'] * test_scaled['measurement_8']\ntest_scaled['feature2'] = test_scaled['measurement_17'] + test_scaled['measurement_6']\ntest_scaled['feature3'] = test_scaled['measurement_10'] / test_scaled['measurement_12']\ntest_scaled['feature4'] = test_scaled['measurement_17'] * test_scaled['measurement_9']","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:51.849473Z","iopub.execute_input":"2022-08-07T00:02:51.849898Z","iopub.status.idle":"2022-08-07T00:02:51.864538Z","shell.execute_reply.started":"2022-08-07T00:02:51.849862Z","shell.execute_reply":"2022-08-07T00:02:51.862806Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# <b>4 <span style='color:#3f4d63'>|</span> Modeling</b>\n\n<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>4.1 | Machine Learning Models</b></p>\n</div>\n\nTo start with, we're gonna use **PyCaret library** in order to examine which algorithms are the ones that we should take into consideration for the submission. In other words, the ones getting the highest **AUC Scores**. ","metadata":{}},{"cell_type":"code","source":"# Reference: https://pycaret.gitbook.io/docs/get-started/functions/train#compare_models\n!pip install --ignore-installed --pre pycaret\nclear_output()\nfrom pycaret.classification import *\nprint('PyCaret setup: ')\nexp_name = setup(data = train_scaled,  target = 'failure')\nprint('\\nComparing models: ')\n# We're going to return top3 models\nbest_model = compare_models(n_select = 3, sort = 'AUC', fold = 5)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:02:51.867080Z","iopub.execute_input":"2022-08-07T00:02:51.868391Z","iopub.status.idle":"2022-08-07T00:08:15.593661Z","shell.execute_reply.started":"2022-08-07T00:02:51.868326Z","shell.execute_reply":"2022-08-07T00:08:15.591949Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Once we have determined which algorithms look more promising, we're gonna tune them. To do so, we can use several tuning libraries. In my case, I'm gonna go for **Optuna**. \n\n#### **Optuna**","metadata":{}},{"cell_type":"code","source":"# Based on MI and Permutation Importance\nlr_features = ['loading','measurement_17','measurement_0','measurement_7','measurement_13','measurement_12',\n              'measurement_3','measurement_8','feature3','feature2','feature4']","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:08:15.595795Z","iopub.execute_input":"2022-08-07T00:08:15.597059Z","iopub.status.idle":"2022-08-07T00:08:15.605954Z","shell.execute_reply.started":"2022-08-07T00:08:15.596991Z","shell.execute_reply":"2022-08-07T00:08:15.603895Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.discriminant_analysis import LinearDiscriminantAnalysis as lda\nfrom sklearn.linear_model import LogisticRegression as lr\nfrom sklearn.ensemble import GradientBoostingClassifier as gbc\nfrom sklearn.metrics import roc_auc_score\nfrom sklearn.model_selection import train_test_split, StratifiedGroupKFold, GroupKFold\n!pip install optuna\nclear_output()\nimport optuna\nfrom optuna.samplers import TPESampler\n\ndef objective(trial):\n    val_scores = []\n    train_scores = []\n    \n    params = {\n        \"penalty\":trial.suggest_categorical(\"penalty\", ['l2']),\n        \"solver\":trial.suggest_categorical(\"solver\", ['liblinear', 'lbfgs', 'newton-cg']),\n        #\"fit_intercept\":trial.suggest_categorical(\"fit_intercept\", [True,False]),\n        \"C\":trial.suggest_float(\"C\",1,50),\n        \"tol\":trial.suggest_float(\"tol\", 1e-5, 1e-2),\n    }\n\n    model = lr(**params)\n    kf = GroupKFold(n_splits=5)   # group by product_code\n\n    for fold, (idx_train, idx_valid) in enumerate(kf.split(train_scaled.drop('failure', axis = 1), \n                                                     train_scaled.failure, groups=train['product_code'])):\n\n        X_train = train_scaled.drop('failure', axis = 1).iloc[idx_train]\n        X_valid = train_scaled.drop('failure', axis = 1).iloc[idx_valid]\n        y_train = train_scaled.loc[idx_train,'failure']\n        y_valid = train_scaled.loc[idx_valid,'failure']\n\n        model.fit(X_train[lr_features], y_train,)        \n        \n        y_train_pred = (model.predict_proba(X_train[lr_features]))[:,1]\n        y_valid_pred = (model.predict_proba(X_valid[lr_features]))[:,1]        \n        train_auc_score = roc_auc_score(y_train, y_train_pred)\n        valid_auc_score = roc_auc_score(y_valid, y_valid_pred)                    \n        train_scores.append(train_auc_score)\n        val_scores.append(valid_auc_score)\n        \n    print(f'Average Accuracy Score of Train: {np.mean(train_scores)}')\n    print(f'Average Accuracy Score of Validation: {np.mean(val_scores)}')\n    return np.mean(val_scores) \n\nallow_optimize = 1\nTRIALS = 150\nTIMEOUT = 3600\n\nif allow_optimize:\n    sampler = TPESampler(seed=123)\n\n    study = optuna.create_study(\n        study_name = 'cat_parameter_opt',\n        direction = 'maximize',\n        sampler = sampler,\n    )\n    study.optimize(objective, n_trials=TRIALS)\n    print(\"Best Score:\",study.best_value)\n    print(\"Best trial\",study.best_trial.params)\n    \n    best_params = study.best_params\n    model_lr = lr(**best_params)\n    \n    kf = StratifiedGroupKFold(n_splits=5)   # group by product_code\n    for fold, (idx_train, idx_valid) in enumerate(kf.split(train_scaled.drop('failure', axis = 1), \n                                                     train_scaled.failure, groups=train['product_code'])):\n\n        X_train = train_scaled.drop('failure', axis = 1).iloc[idx_train]\n        X_valid = train_scaled.drop('failure', axis = 1).iloc[idx_valid]\n        y_train = train_scaled.loc[idx_train,'failure']\n        y_valid = train_scaled.loc[idx_valid,'failure']\n\n        model_lr.fit(X_train[lr_features], y_train)\nelse:\n    X_train_tmp, X_valid_tmp, y_train_tmp, y_valid_tmp = train_test_split(X, y, test_size=0.3, random_state=42)\n    model_cat = CatBoostClassifier(\n        verbose=1000,\n        early_stopping_rounds=10,\n        #iterations=5000,\n        random_state = 2022, learning_rate = 0.08665686887824392, bagging_temperature = 2.010272294890727, n_estimators = 806, max_depth = 7, \n        random_strength = 35, l2_leaf_reg = 1.2373460332766636e-05, min_child_samples = 69, max_bin = 317, od_type = 'IncToDec', \n        task_type = 'GPU', eval_metric = 'Accuracy'\n    ).fit(X_train_tmp, y_train_tmp, eval_set=[(X_valid_tmp, y_valid_tmp)], early_stopping_rounds=35)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:08:15.608953Z","iopub.execute_input":"2022-08-07T00:08:15.609646Z","iopub.status.idle":"2022-08-07T00:10:12.067693Z","shell.execute_reply.started":"2022-08-07T00:08:15.609449Z","shell.execute_reply":"2022-08-07T00:10:12.066395Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_feature_importance(model_lr.coef_.flatten(), lr_features,'Logistic Regression')    ","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:10:36.121388Z","iopub.execute_input":"2022-08-07T00:10:36.121846Z","iopub.status.idle":"2022-08-07T00:10:36.199245Z","shell.execute_reply.started":"2022-08-07T00:10:36.121811Z","shell.execute_reply":"2022-08-07T00:10:36.198099Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import eli5\nfrom eli5.sklearn import PermutationImportance\n\nperm = PermutationImportance(model_lr, random_state=1).fit(train_scaled.drop('failure', axis = 1)[lr_features], \n                                                     train_scaled.failure)\neli5.show_weights(perm, feature_names = lr_features)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:10:52.031654Z","iopub.execute_input":"2022-08-07T00:10:52.032424Z","iopub.status.idle":"2022-08-07T00:10:58.871303Z","shell.execute_reply.started":"2022-08-07T00:10:52.032370Z","shell.execute_reply":"2022-08-07T00:10:58.869639Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"preds = (model_lr.predict_proba(test_scaled[lr_features]))[:,1]\nss_lr = sample_submission.copy()\nss_lr['failure'] = preds\n#ss_lr.to_csv('submission.csv', index=False)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:11:20.090796Z","iopub.execute_input":"2022-08-07T00:11:20.091210Z","iopub.status.idle":"2022-08-07T00:11:20.106942Z","shell.execute_reply.started":"2022-08-07T00:11:20.091178Z","shell.execute_reply":"2022-08-07T00:11:20.105270Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def objective(trial):\n    val_scores = []\n    train_scores = []\n    \n    params = {\n        \"tol\":trial.suggest_float(\"tol\", 1e-5, 1e-2),\n    }\n\n    model = lda(**params)\n    kf = StratifiedGroupKFold(n_splits=5)   # group by product_code\n\n    for fold, (idx_train, idx_valid) in enumerate(kf.split(train_scaled.drop('failure', axis = 1), \n                                                     train_scaled.failure, groups=train['product_code'])):\n\n        X_train = train_scaled.drop('failure', axis = 1).iloc[idx_train]\n        X_valid = train_scaled.drop('failure', axis = 1).iloc[idx_valid]\n        y_train = train_scaled.loc[idx_train,'failure']\n        y_valid = train_scaled.loc[idx_valid,'failure']\n\n        model.fit(X_train, y_train,)        \n        \n        y_train_pred = (model.predict_proba(X_train))[:,1]\n        y_valid_pred = (model.predict_proba(X_valid))[:,1]\n        train_auc_score = roc_auc_score(y_train, y_train_pred)\n        valid_auc_score = roc_auc_score(y_valid, y_valid_pred)                    \n        train_scores.append(train_auc_score)\n        val_scores.append(valid_auc_score)\n        \n    print(f'Average Accuracy Score of Train: {np.mean(train_scores)}')\n    print(f'Average Accuracy Score of Validation: {np.mean(val_scores)}')\n    return np.mean(val_scores)\n\nallow_optimize = 1\nTRIALS = 150\nTIMEOUT = 3600\n\nif allow_optimize:\n    sampler = TPESampler(seed=123)\n\n    study = optuna.create_study(\n        study_name = 'cat_parameter_opt',\n        direction = 'maximize',\n        sampler = sampler,\n    )\n    study.optimize(objective, n_trials=TRIALS)\n    print(\"Best Score:\",study.best_value)\n    print(\"Best trial\",study.best_trial.params)\n    \n    best_params = study.best_params\n    model_lda = lda(**best_params)\n    \n    kf = StratifiedGroupKFold(n_splits=5)   # group by product_code\n    for fold, (idx_train, idx_valid) in enumerate(kf.split(train_scaled.drop('failure', axis = 1), \n                                                     train_scaled.failure, groups=train['product_code'])):\n\n        X_train = train_scaled.drop('failure', axis = 1).iloc[idx_train]\n        X_valid = train_scaled.drop('failure', axis = 1).iloc[idx_valid]\n        y_train = train_scaled.loc[idx_train,'failure']\n        y_valid = train_scaled.loc[idx_valid,'failure']\n\n        model_lda.fit(X_train, y_train)","metadata":{"_kg_hide-output":true,"execution":{"iopub.status.busy":"2022-08-07T00:11:24.045495Z","iopub.execute_input":"2022-08-07T00:11:24.046732Z","iopub.status.idle":"2022-08-07T00:14:04.747729Z","shell.execute_reply.started":"2022-08-07T00:11:24.046676Z","shell.execute_reply":"2022-08-07T00:14:04.746460Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_feature_importance(model_lda.coef_.flatten(), train_scaled.drop('failure',axis=1).columns,'Logistic Regression')    ","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:14:04.750493Z","iopub.execute_input":"2022-08-07T00:14:04.758130Z","iopub.status.idle":"2022-08-07T00:14:04.894486Z","shell.execute_reply.started":"2022-08-07T00:14:04.758008Z","shell.execute_reply":"2022-08-07T00:14:04.893243Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"preds = (model_lda.predict_proba(test_scaled))[:,1]\nss_lda = sample_submission.copy()\nss_lda['failure'] = preds","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:14:04.896392Z","iopub.execute_input":"2022-08-07T00:14:04.898401Z","iopub.status.idle":"2022-08-07T00:14:04.912637Z","shell.execute_reply.started":"2022-08-07T00:14:04.898361Z","shell.execute_reply":"2022-08-07T00:14:04.910773Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<div style=\"color:white;display:fill;\n            background-color:#3f4d6f;font-size:150%;\n            font-family:Nexa;letter-spacing:0.5px\">\n    <p style=\"padding: 8px;color:white;\"><b>4.2 | Neural Network. Keras approach</b></p>\n</div>","metadata":{}},{"cell_type":"code","source":"import tensorflow as tf\nfrom tensorflow import keras\nfrom tensorflow.keras import layers\nfrom tensorflow.keras.utils import plot_model\nfrom tensorflow.keras.callbacks import EarlyStopping\n\ninitializer = tf.keras.initializers.HeNormal()\nmodel = keras.Sequential([\n    # ====== Input Layer ======\n    layers.Input(shape=[train_scaled.drop('failure',axis=1).shape[1]]),\n    \n    # ======  Hidden ReLU layers ======\n    layers.Dropout(rate=0.3), # apply 30% dropout to the next layer\n    layers.BatchNormalization(),\n    layers.Dense(units=128, activation='relu', kernel_initializer = initializer),\n    \n    layers.Dropout(rate=0.2), \n    layers.BatchNormalization(),\n    layers.Dense(units=32, activation='swish', kernel_initializer = initializer),\n    \n    layers.Dropout(rate=0.15), \n    layers.BatchNormalization(),    \n    layers.Dense(units=32, activation='swish', kernel_initializer = initializer),\n    \n    layers.Dropout(rate=0.05), \n    layers.BatchNormalization(),    \n    layers.Dense(units=8, activation='swish', kernel_initializer = initializer),\n    \n    # ====== Output layer ====== \n    # Binary Classification -> sigmoid activation function\n    layers.Dense(units=1, activation = 'sigmoid'),\n])\n\nplot_model(model, show_shapes=True)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-07T00:14:04.923166Z","iopub.execute_input":"2022-08-07T00:14:04.928355Z","iopub.status.idle":"2022-08-07T00:14:07.865655Z","shell.execute_reply.started":"2022-08-07T00:14:04.928270Z","shell.execute_reply":"2022-08-07T00:14:07.864319Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.model_selection import train_test_split, StratifiedGroupKFold\n\nLOSS = tf.keras.losses.BinaryCrossentropy()\nMETRIC = tf.keras.metrics.AUC(name='auc')\noptimizer = tf.keras.optimizers.Adam(learning_rate = 0.00125)\nmodel.compile(\n    optimizer=optimizer,\n    loss=LOSS,\n    metrics=METRIC,\n)\n\nearly_stopping = EarlyStopping(\n    min_delta=0.0001, # minimium amount of change to count as an improvement\n    patience=20, # how many epochs to wait before stopping\n    restore_best_weights=True,\n)\n\n# ================================\n# Code taken from: https://www.kaggle.com/code/pourchot/keras-optuna-for-logistic-reg\n# Author: @pourchot\ncheckpoint_filepath = 'model'\nCheckpoint = tf.keras.callbacks.ModelCheckpoint(\n    filepath = checkpoint_filepath,\n    monitor=\"val_auc\",\n    verbose=0,\n    save_best_only=True,\n    save_weights_only=True,\n    mode=\"max\",\n    save_freq=\"epoch\",\n    options=None,\n    initial_value_threshold=None\n)\n# ================================\n\nkf = GroupKFold(n_splits=5)   # group by product_code\n\nfor fold, (idx_train, idx_valid) in enumerate(kf.split(train_scaled.drop('failure', axis = 1), \n                                                 train_scaled.failure, groups=train_scaled['product_code'])):\n    \n    X_train = train_scaled.drop('failure', axis = 1).iloc[idx_train]\n    X_valid = train_scaled.drop('failure', axis = 1).iloc[idx_valid]\n    y_train = train_scaled.loc[idx_train,'failure']\n    y_valid = train_scaled.loc[idx_valid,'failure']\n\n    history = model.fit(\n            X_train, y_train,\n            validation_data=(X_valid, y_valid),\n            batch_size=128,\n            epochs=100,\n            callbacks=[early_stopping, Checkpoint],\n            #verbose=0, # hide the output because we have so many epochs\n        )        \n    \n    history_df = pd.DataFrame(history.history)\n    history_df.loc[:, ['loss','val_loss']].plot()\n    plt.ylim(0.45, 0.65) \n    \n    history_df.loc[:, ['auc','val_auc']].plot()\n    plt.ylim(0.475, 0.725)        \n    # We plot the epoch in which we achive maximum rocauc score\n    max_idx = np.argmax(history_df['val_auc'])\n    plt.scatter(x = max_idx, y = history_df.loc[max_idx, 'val_auc'], color = 'red')","metadata":{"_kg_hide-output":true,"execution":{"iopub.status.busy":"2022-08-07T00:14:07.867740Z","iopub.execute_input":"2022-08-07T00:14:07.869363Z","iopub.status.idle":"2022-08-07T00:17:04.331961Z","shell.execute_reply.started":"2022-08-07T00:14:07.869310Z","shell.execute_reply":"2022-08-07T00:17:04.330278Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"preds = model.predict(test_scaled)\nss_nn = sample_submission.copy()\nss_nn['failure'] = preds\n#ss_nn.to_csv('submission.csv', index=False)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:17:04.333965Z","iopub.execute_input":"2022-08-07T00:17:04.334453Z","iopub.status.idle":"2022-08-07T00:17:05.651090Z","shell.execute_reply.started":"2022-08-07T00:17:04.334415Z","shell.execute_reply":"2022-08-07T00:17:05.650112Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.distplot(ss_lr['failure'], label = 'Logistic Regression')\nsns.distplot(ss_lda['failure'], label = 'LDA')\nsns.distplot(ss_nn['failure'], label = 'Neural Network')\nplt.legend()","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:17:05.652301Z","iopub.execute_input":"2022-08-07T00:17:05.652802Z","iopub.status.idle":"2022-08-07T00:17:06.561775Z","shell.execute_reply.started":"2022-08-07T00:17:05.652770Z","shell.execute_reply":"2022-08-07T00:17:06.560594Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ss = sample_submission.copy()\nss['failure'] = 1 * ss_lr['failure'] + ss_nn['failure'] * 0\nss.to_csv('submission.csv', index=False)","metadata":{"execution":{"iopub.status.busy":"2022-08-07T00:17:06.563377Z","iopub.execute_input":"2022-08-07T00:17:06.564573Z","iopub.status.idle":"2022-08-07T00:17:06.629398Z","shell.execute_reply.started":"2022-08-07T00:17:06.564534Z","shell.execute_reply":"2022-08-07T00:17:06.627721Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## References\n\n* [Keras + Optuna for Logistic Regression. Author: @pourchot](https://www.kaggle.com/code/pourchot/keras-optuna-for-logistic-reg)\n* [TPS Aug22 - Failure Prediction. Author: @samuelcortinhas](https://www.kaggle.com/code/samuelcortinhas/tps-aug-22-failure-prediction)\n* [TPSAUG22 EDA which makes sense. Author: @ambrosm](https://www.kaggle.com/code/ambrosm/tpsaug22-eda-which-makes-sense)","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}