{"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":"# Missing Value Fill Deployment\n\nupdates from the very interesting notebook : https://www.kaggle.com/code/takanashihumbert/tps-aug22-lb-0-59013","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\n\n\nfrom matplotlib.ticker import MaxNLocator\nimport seaborn as sns\nfrom cycler import cycler\nfrom IPython.display import display\nimport math\nimport os\nimport random\nimport gc\nimport sys\nimport warnings\nwarnings.filterwarnings('ignore')\n\nimport optuna\nfrom colorama import Fore, Back, Style\n\nfrom sklearn.model_selection import StratifiedGroupKFold, StratifiedKFold, train_test_split,GroupKFold\nfrom sklearn.metrics import roc_auc_score\nfrom sklearn.calibration import CalibrationDisplay\nfrom sklearn.preprocessing import StandardScaler,RobustScaler,LabelEncoder\nfrom sklearn.impute import KNNImputer\nfrom sklearn import linear_model\nfrom sklearn.linear_model import HuberRegressor\nfrom sklearn.decomposition import PCA\nfrom sklearn.naive_bayes import BernoulliNB\nfrom sklearn.neighbors import KNeighborsClassifier\nfrom sklearn.calibration import CalibratedClassifierCV\n\nimport matplotlib.pyplot as plt\n\nimport tensorflow as tf\nfrom tensorflow import keras\nfrom tensorflow.keras import layers\nfrom tensorflow.keras.callbacks import ReduceLROnPlateau, LearningRateScheduler, EarlyStopping\nfrom tensorflow.keras.layers import Input, Dense, Activation,  BatchNormalization, Dropout, Concatenate, Embedding,  Flatten, Conv1D\nfrom tensorflow.keras.models import Model","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:07.972970Z","iopub.execute_input":"2022-08-12T03:35:07.973391Z","iopub.status.idle":"2022-08-12T03:35:07.984959Z","shell.execute_reply.started":"2022-08-12T03:35:07.973359Z","shell.execute_reply":"2022-08-12T03:35:07.983528Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train = pd.read_csv('../input/tabular-playground-series-aug-2022/train.csv')\ntest = pd.read_csv('../input/tabular-playground-series-aug-2022/test.csv')\nsubmission = pd.read_csv('../input/tabular-playground-series-aug-2022/sample_submission.csv')\ntarget = train['failure']\ntrain.drop('failure',axis=1, inplace = True)\ndata = pd.concat([train, test])\ntrain.shape,test.shape","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:08.281753Z","iopub.execute_input":"2022-08-12T03:35:08.282417Z","iopub.status.idle":"2022-08-12T03:35:08.600675Z","shell.execute_reply.started":"2022-08-12T03:35:08.282382Z","shell.execute_reply":"2022-08-12T03:35:08.599842Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Feature Engineering","metadata":{}},{"cell_type":"markdown","source":"\nInitial preprocessing from :\n\nhttps://www.kaggle.com/code/desalegngeb/tps08-logisticregression-and-some-fe\n\nhttps://www.kaggle.com/code/alikayed/tps08-logisticregression-and-some-fe-c83a47","metadata":{}},{"cell_type":"code","source":"# library for coding string values :\n! pip install feature_engine\nfrom feature_engine.encoding import WoEEncoder","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:08.972109Z","iopub.execute_input":"2022-08-12T03:35:08.972525Z","iopub.status.idle":"2022-08-12T03:35:23.918425Z","shell.execute_reply.started":"2022-08-12T03:35:08.972492Z","shell.execute_reply":"2022-08-12T03:35:23.917132Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data['m3_missing'] = data['measurement_3'].isnull().astype(np.int8)\ndata['m5_missing'] = data['measurement_5'].isnull().astype(np.int8)\ndata['area'] = data['attribute_2'] * data['attribute_3']\n\nfeature = [f for f in test.columns if f.startswith('measurement') or f=='loading']\n\n# dictionnary of dictionnaries (for the 11 best correlated measurement columns), \n# we will use the dictionnaries below to select the best correlated columns according to the product code)\n# Only for 'measurement_17' we make a 'manual' selection :\nfull_fill_dict ={}\nfull_fill_dict['measurement_17'] = {\n    'A': ['measurement_5','measurement_6','measurement_8'],\n    'B': ['measurement_4','measurement_5','measurement_7'],\n    'C': ['measurement_5','measurement_7','measurement_8','measurement_9'],\n    'D': ['measurement_5','measurement_6','measurement_7','measurement_8'],\n    'E': ['measurement_4','measurement_5','measurement_6','measurement_8'],\n    'F': ['measurement_4','measurement_5','measurement_6','measurement_7'],\n    'G': ['measurement_4','measurement_6','measurement_8','measurement_9'],\n    'H': ['measurement_4','measurement_5','measurement_7','measurement_8','measurement_9'],\n    'I': ['measurement_3','measurement_7','measurement_8']\n}\n\n# collect the name of the next 10 best measurement columns sorted by correlation (except 17 already done above):\ncol = [col for col in test.columns if 'measurement' not in col]+ ['loading','m3_missing','m5_missing']\na = []\nb =[]\nfor x in range(3,17):\n    corr = np.absolute(data.drop(col, axis=1).corr()[f'measurement_{x}']).sort_values(ascending=False)\n    a.append(np.round(np.sum(corr[1:4]),3)) # we add the 3 first lines of the correlation values to get the \"most correlated\"\n    b.append(f'measurement_{x}')\nc = pd.DataFrame()\nc['Selected columns'] = b\nc['correlation total'] = a\nc = c.sort_values(by = 'correlation total',ascending=False).reset_index(drop = True)\nprint(f'Columns selected by correlation sum of the 3 first rows : ')\ndisplay(c.head(10))\n\nfor i in range(10):\n    measurement_col = 'measurement_' + c.iloc[i,0][12:] # we select the next best correlated column \n    fill_dict ={}\n    for x in data.product_code.unique() : \n        corr = np.absolute(data[data.product_code == x].drop(col, axis=1).corr()[measurement_col]).sort_values(ascending=False)\n        measurement_col_dic = {}\n        measurement_col_dic[measurement_col] = corr[1:5].index.tolist()\n        fill_dict[x] = measurement_col_dic[measurement_col]\n    full_fill_dict[measurement_col] =fill_dict\n    \nfeature = [f for f in data.columns if f.startswith('measurement') or f=='loading']\nnullValue_cols = [col for col in train.columns if train[col].isnull().sum()!=0]\n    \nfor code in data.product_code.unique():\n    total_na_filled_by_linear_model = 0\n    print(f'\\n-------- Product code {code} ----------\\n')\n    print(f'filled by linear model :')\n    for measurement_col in list(full_fill_dict.keys()):\n        tmp = data[data.product_code==code]\n        column = full_fill_dict[measurement_col][code]\n        tmp_train = tmp[column+[measurement_col]].dropna(how='any')\n        tmp_test = tmp[(tmp[column].isnull().sum(axis=1)==0)&(tmp[measurement_col].isnull())]\n\n        model = HuberRegressor(epsilon=1.9)\n        model.fit(tmp_train[column], tmp_train[measurement_col])\n        data.loc[(data.product_code==code)&(data[column].isnull().sum(axis=1)==0)&(data[measurement_col].isnull()),measurement_col] = model.predict(tmp_test[column])\n        print(f'{measurement_col} : {len(tmp_test)}')\n        total_na_filled_by_linear_model += len(tmp_test)\n        \n    # others NA columns:\n    NA = data.loc[data[\"product_code\"] == code,nullValue_cols ].isnull().sum().sum()\n    model1 = KNNImputer(n_neighbors=3)\n    data.loc[data.product_code==code, feature] = model1.fit_transform(data.loc[data.product_code==code, feature])\n    print(f'\\n{total_na_filled_by_linear_model} filled by linear model ') \n    print(f'{NA} filled by KNN ')\n    \ndata['measurement_avg'] = data[[f'measurement_{i}' for i in range(3, 17)]].mean(axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:23.921396Z","iopub.execute_input":"2022-08-12T03:35:23.921994Z","iopub.status.idle":"2022-08-12T03:35:46.192386Z","shell.execute_reply.started":"2022-08-12T03:35:23.921944Z","shell.execute_reply":"2022-08-12T03:35:46.191181Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def _scale(train_data, val_data, test_data, feats):\n    scaler = StandardScaler()\n    # scaler = PowerTransformer()\n    \n    scaled_train = scaler.fit_transform(train_data[feats])\n    scaled_val = scaler.transform(val_data[feats])\n    scaled_test = scaler.transform(test_data[feats])\n    \n    #back to dataframe\n    new_train = train_data.copy()\n    new_val = val_data.copy()\n    new_test = test_data.copy()\n    \n    new_train[feats] = scaled_train\n    new_val[feats] = scaled_val\n    new_test[feats] = scaled_test\n    \n    assert len(train_data) == len(new_train)\n    assert len(val_data) == len(new_val)\n    assert len(test_data) == len(new_test)\n    \n    return new_train, new_val, new_test","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:46.193970Z","iopub.execute_input":"2022-08-12T03:35:46.194422Z","iopub.status.idle":"2022-08-12T03:35:46.202957Z","shell.execute_reply.started":"2022-08-12T03:35:46.194389Z","shell.execute_reply":"2022-08-12T03:35:46.201535Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train = data.iloc[:train.shape[0],:]\ntest = data.iloc[train.shape[0]:,:]\nprint(train.shape, test.shape)\n\ngroups = train.product_code\nX = train\ny = target","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:46.206211Z","iopub.execute_input":"2022-08-12T03:35:46.206597Z","iopub.status.idle":"2022-08-12T03:35:46.227182Z","shell.execute_reply.started":"2022-08-12T03:35:46.206563Z","shell.execute_reply":"2022-08-12T03:35:46.226137Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Thanks to @MAXSARMENTO \nwoe_encoder = WoEEncoder(variables=['attribute_0'])\nwoe_encoder.fit(X, y)\nX = woe_encoder.transform(X)\ntest = woe_encoder.transform(test)","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:46.228298Z","iopub.execute_input":"2022-08-12T03:35:46.228785Z","iopub.status.idle":"2022-08-12T03:35:46.306360Z","shell.execute_reply.started":"2022-08-12T03:35:46.228715Z","shell.execute_reply":"2022-08-12T03:35:46.305446Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"select_feature = ['loading',\n                  'attribute_0',\n                  'measurement_17',\n                  'measurement_0',\n                  'measurement_1',\n                  'measurement_2',\n                  'area',\n                  'm3_missing',\n                  'm5_missing',\n                  'measurement_avg']","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:46.308134Z","iopub.execute_input":"2022-08-12T03:35:46.308958Z","iopub.status.idle":"2022-08-12T03:35:46.314640Z","shell.execute_reply.started":"2022-08-12T03:35:46.308913Z","shell.execute_reply":"2022-08-12T03:35:46.313784Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## VAE Implementation","metadata":{}},{"cell_type":"code","source":"from keras.layers import Lambda, Input, Dense\nfrom keras.models import Model\nfrom keras.datasets import mnist\nfrom keras.losses import mse, binary_crossentropy\nfrom keras import backend as K\nfrom keras.callbacks import ModelCheckpoint\nfrom keras.layers import Input, Dense, Lambda, Layer, Add, Multiply\nfrom keras.models import Model, Sequential\nfrom sklearn.preprocessing import MinMaxScaler","metadata":{"execution":{"iopub.status.busy":"2022-08-12T04:21:21.012953Z","iopub.execute_input":"2022-08-12T04:21:21.013418Z","iopub.status.idle":"2022-08-12T04:21:21.020728Z","shell.execute_reply.started":"2022-08-12T04:21:21.013380Z","shell.execute_reply":"2022-08-12T04:21:21.019283Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"original_dim= len(select_feature)\ninput_shape = (original_dim, )\nintermediate_dim = int(original_dim/2)\nbatch_size = 256\nlatent_dim = 10\nepochs     = 300\nepsilon_std = 1.0","metadata":{"execution":{"iopub.status.busy":"2022-08-12T04:21:25.728806Z","iopub.execute_input":"2022-08-12T04:21:25.729219Z","iopub.status.idle":"2022-08-12T04:21:25.735105Z","shell.execute_reply.started":"2022-08-12T04:21:25.729187Z","shell.execute_reply":"2022-08-12T04:21:25.734005Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class KLDivergenceLayer(Layer):\n\n    \"\"\" Identity transform layer that adds KL divergence\n    to the final model loss.\n    \"\"\"\n\n    def __init__(self, *args, **kwargs):\n        self.is_placeholder = True\n        super(KLDivergenceLayer, self).__init__(*args, **kwargs)\n\n    def call(self, inputs):\n\n        mu, log_var = inputs\n\n        kl_batch = - .5 * K.sum(1 + log_var -\n                                K.square(mu) -\n                                K.exp(log_var), axis=-1)\n\n        self.add_loss(K.mean(kl_batch), inputs=inputs)\n\n        return inputs","metadata":{"execution":{"iopub.status.busy":"2022-08-12T04:21:26.673006Z","iopub.execute_input":"2022-08-12T04:21:26.673449Z","iopub.status.idle":"2022-08-12T04:21:26.680935Z","shell.execute_reply.started":"2022-08-12T04:21:26.673414Z","shell.execute_reply":"2022-08-12T04:21:26.679838Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# VAE Architecture\n# * original_dim - Original Input Dimension\n# * intermediate_dim - Hidden Layer Dimension\n# * latent_dim - Latent/Embedding Dimension\ndef vae_arc(original_dim, intermediate_dim, latent_dim):\n    # Decode\n    decoder = Sequential([\n        Dense(intermediate_dim, input_dim=latent_dim, activation='relu'),\n        Dense(original_dim, activation='sigmoid')\n    ])\n\n    # Encode\n    x = Input(shape=(original_dim,))\n    h = Dense(intermediate_dim, activation='relu')(x)\n\n    z_mu = Dense(latent_dim)(h)\n    z_log_var = Dense(latent_dim)(h)\n\n    z_mu, z_log_var = KLDivergenceLayer()([z_mu, z_log_var])\n    z_sigma = Lambda(lambda t: K.exp(.5*t))(z_log_var)\n\n    eps = Input(tensor=K.random_normal(stddev=epsilon_std,\n                                       shape=(K.shape(x)[0], latent_dim)))\n    z_eps = Multiply()([z_sigma, eps])\n    z = Add()([z_mu, z_eps])\n\n    x_pred = decoder(z)\n    \n    return x, eps, z_mu, x_pred","metadata":{"execution":{"iopub.status.busy":"2022-08-12T04:21:28.131494Z","iopub.execute_input":"2022-08-12T04:21:28.132722Z","iopub.status.idle":"2022-08-12T04:21:28.141627Z","shell.execute_reply.started":"2022-08-12T04:21:28.132683Z","shell.execute_reply":"2022-08-12T04:21:28.140414Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def nll(y_true, y_pred):\n    \"\"\" Negative log likelihood (Bernoulli). \"\"\"\n\n    # keras.losses.binary_crossentropy gives the mean\n    # over the last axis. we require the sum\n    return K.sum(K.binary_crossentropy(y_true, y_pred), axis=-1)","metadata":{"execution":{"iopub.status.busy":"2022-08-12T04:21:28.911143Z","iopub.execute_input":"2022-08-12T04:21:28.912032Z","iopub.status.idle":"2022-08-12T04:21:28.918293Z","shell.execute_reply.started":"2022-08-12T04:21:28.911976Z","shell.execute_reply":"2022-08-12T04:21:28.917164Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"x, eps, z_mu, x_pred = vae_arc(original_dim, intermediate_dim, latent_dim)","metadata":{"execution":{"iopub.status.busy":"2022-08-12T04:25:53.069623Z","iopub.execute_input":"2022-08-12T04:25:53.070091Z","iopub.status.idle":"2022-08-12T04:25:53.163709Z","shell.execute_reply.started":"2022-08-12T04:25:53.070046Z","shell.execute_reply":"2022-08-12T04:25:53.162447Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"callbacks_list = [tf.keras.callbacks.EarlyStopping(\n    monitor=\"val_loss\",\n    min_delta=0,\n    patience=10,\n    verbose=0,\n    mode=\"auto\",\n    baseline=None,\n    restore_best_weights=True,\n)]","metadata":{"execution":{"iopub.status.busy":"2022-08-12T04:21:34.381643Z","iopub.execute_input":"2022-08-12T04:21:34.382876Z","iopub.status.idle":"2022-08-12T04:21:34.389424Z","shell.execute_reply.started":"2022-08-12T04:21:34.382818Z","shell.execute_reply":"2022-08-12T04:21:34.388191Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Training","metadata":{}},{"cell_type":"code","source":"print(f'\\n******** cross validation strategy : Without folds ********\\n')\n\nlr_test = np.zeros(len(test))\nimportance_list = []\n\nmodel = linear_model.LogisticRegression(max_iter=200, C=0.0001, penalty='l2', solver='newton-cg')\nmodel.fit(X[select_feature], y)\nimportance_list.append(model.coef_.ravel())\n\nlr_test += model.predict_proba(test[select_feature])[:, 1]\n\nimportance_df = pd.DataFrame(np.array(importance_list).T, index=X[select_feature].columns)\nimportance_df['mean'] = importance_df.mean(axis=1).abs()\nimportance_df['feature'] = X[select_feature].columns\nimportance_df = importance_df.sort_values('mean', ascending=False).reset_index().head(20)\nplt.barh(importance_df.index, importance_df['mean'], color='lightgreen')\nplt.gca().invert_yaxis()\nplt.yticks(ticks=importance_df.index, labels=importance_df['feature'])\nplt.title('LogisticRegression feature importances')\nplt.show()\n\nsubmission['failure'] = lr_test\nsubmission.to_csv(f\"./sub_missing_no_fold.csv\", index=False)\nsubmission\n\nfor step,cross_val in enumerate([StratifiedKFold(n_splits=5, shuffle=True, random_state=0),GroupKFold(n_splits=5)]) :\n    print(f'\\n******** cross validation strategy : {cross_val} ********\\n')\n    lr_oof = np.zeros(len(train))\n    lr_test = np.zeros(len(test))\n    lr_auc = 0\n    importance_list = []\n    \n    kf = cross_val\n    for fold_idx, (train_idx, val_idx) in enumerate(kf.split(X, y,groups = train.product_code )):\n        x_train, x_val = X.iloc[train_idx], X.iloc[val_idx]\n        y_train, y_val = y.iloc[train_idx], y.iloc[val_idx]\n        x_train, x_val, x_test = _scale(x_train, x_val, test, select_feature)\n        \n        scaler    = MinMaxScaler()\n        X_train_norm   = scaler.fit_transform(x_train[select_feature])\n        X_val_norm = scaler.transform(x_val[select_feature])\n        X_test_norm = scaler.transform(x_test[select_feature])\n        \n        vae = Model(inputs=[x, eps], outputs=x_pred)\n        vae.compile(optimizer='adam', loss=nll)\n        \n        vae.fit([X_train_norm, X_train_norm],[X_train_norm, X_train_norm],\n                        epochs=epochs,\n                        batch_size=batch_size,\n                        verbose=0,\n                        callbacks=callbacks_list)\n        \n        X_train_vae = vae.predict([X_train_norm, X_train_norm])\n        X_test_vae = vae.predict([X_test_norm, X_test_norm])\n        X_val_vae = vae.predict([X_val_norm, X_val_norm])\n\n        model = linear_model.LogisticRegression(max_iter=200, C=0.0001, penalty='l2', solver='newton-cg')\n        model.fit(X_train_vae, y_train)\n        importance_list.append(model.coef_.ravel())\n\n        val_preds = model.predict_proba(X_val_vae)[:, 1]\n        print(\"FOLD: \", fold_idx+1, \" ROC-AUC:\", round(roc_auc_score(y_val, val_preds), 5))\n        lr_auc += roc_auc_score(y_val, val_preds) / 5\n        lr_test += model.predict_proba(X_test_vae)[:, 1] / 5\n        lr_oof[val_idx] = val_preds\n\n    print(f\"\\n{Fore.GREEN}{Style.BRIGHT}Average auc = {round(lr_auc, 5)}{Style.RESET_ALL}\")\n    print(f\"{Fore.BLUE}{Style.BRIGHT}OOF auc     = {round(roc_auc_score(y, lr_oof), 5)}{Style.RESET_ALL}\\n\")\n\n    importance_df = pd.DataFrame(np.array(importance_list).T, index=x_train[select_feature].columns)\n    importance_df['mean'] = importance_df.mean(axis=1).abs()\n    importance_df['feature'] = x_train[select_feature].columns\n    importance_df = importance_df.sort_values('mean', ascending=False).reset_index().head(20)\n    plt.barh(importance_df.index, importance_df['mean'], color='lightgreen')\n    plt.gca().invert_yaxis()\n    plt.yticks(ticks=importance_df.index, labels=importance_df['feature'])\n    plt.title('LogisticRegression feature importances')\n    plt.show()\n \n    submission['failure'] = lr_test\n    submission.to_csv(f\"./sub_missing_{step}.csv\", index=False)\n    submission","metadata":{"execution":{"iopub.status.busy":"2022-08-12T03:35:58.008016Z","iopub.execute_input":"2022-08-12T03:35:58.008460Z","iopub.status.idle":"2022-08-12T03:36:02.648752Z","shell.execute_reply.started":"2022-08-12T03:35:58.008429Z","shell.execute_reply":"2022-08-12T03:36:02.647463Z"},"trusted":true},"execution_count":null,"outputs":[]}]}