{"cells":[{"metadata":{},"cell_type":"markdown","source":"## Classifier Accuracy Across TTF Range\n---\n\nThe purpose of this simple kernel is to evaluate how accurately a series of LightGBM classifiers can predict that `time_to_failure` lies within a specified range, and how this accuracy is affected by the size of TTF. Which are the easier and harder parts of the data to successfully predict? More importantly, I wanted to investigate a completely new approach to this competition. \n\nSo far, every kernel I have seen has worked on the principal of minimising MAE within their ML model, in accordance with the competition objectives. Is it possible that by reframing the nature of our target, we can indirectly achieve better results?\n\nFor simplicity and reproducibility, this kernel uses Andrew's Features only. You may find better results with your own data."},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"import numpy as np \nimport pandas as pd \nimport lightgbm as lgb\nfrom sklearn.model_selection import KFold\nfrom sklearn.metrics import roc_auc_score, mean_absolute_error\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport warnings\nwarnings.filterwarnings('ignore')","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":true},"cell_type":"code","source":"#I took these from a private kernel of mine but you can add Andrew's Features yourself\ntrain = pd.read_csv('../input/andrews-features/train_X.csv')\ntest = pd.read_csv('../input/andrews-features/test_X.csv')\ny_train = pd.read_csv('../input/andrews-features/train_y.csv').values.flatten()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"This quick section identifies the 'overlaps' between earthquake periods so these segments can be removed. It also creates a sequential column for the earthquake number that you can use for subsetting."},{"metadata":{"trusted":true},"cell_type":"code","source":"#remove overlapping segments\ndef sequential(df, col='time_to_failure', newcol='quake_id', overlap='remove_idx'):\n    df[newcol] = np.zeros(len(df))\n    df[overlap] = np.zeros(len(df))\n    for i in range(1, len(df)):\n        if df.loc[i, col] > df.loc[i-1, col]:\n            df.loc[i, newcol] = df.loc[i-1, newcol] + 1\n            df.loc[i, overlap] = 1\n        else:\n            df.loc[i, newcol] = df.loc[i-1, newcol]\n    return(df)   \n\ntrain['time_to_failure'] = y_train\ntrain = sequential(train)\nprint(train.quake_id.describe())\nprint('Total number of overlapping segments: ', train.remove_idx.sum())","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"17 distinct quakes and 16 overlapping segments, as expected. The overlapping rows can be removed, and the array of quake ids kept. You can repeat this process for your own work or just copy the subsetting arrays from this kernel."},{"metadata":{"trusted":true},"cell_type":"code","source":"#keep only non-overlapping segments\nkeep_index = train.loc[train.remove_idx != 1, :].index\ntrain = train.iloc[keep_index].reset_index(drop=True)\ny_train = y_train[keep_index]\n\n#save quake ids as numpy array and remove unnecessary columns\nquake_ids = train['quake_id'].values\nnp.save('quake_ids.npy', quake_ids)\nnp.save('keep_index.npy', keep_index)\ntrain.drop(['remove_idx', 'time_to_failure', 'quake_id'], axis=1, inplace=True)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Now we can run a 5-fold classifier to identify whether TTF lies in a given range of 0.5 seconds. We'll run this function in a loop so the probability of TTF is evaluated for each row across the entire TTF range in train.csv "},{"metadata":{"trusted":true},"cell_type":"code","source":"N_FOLD = 5\nSEP = 0.5\n\nfolds = KFold(n_splits=N_FOLD, shuffle=True, random_state=42)\nfeature_importance = np.zeros(len(train.columns))\n\ntrain_preds = pd.DataFrame()\ntest_preds = pd.DataFrame()\n\ndef ttf_classifier(threshold, X, X_test, df, test_df, y=y_train, sep=SEP, feature_importance=feature_importance):\n    \n    #y == 1 if TTF lies in specific range\n    y = np.logical_and(y >= threshold, y < threshold + sep)\n    models = []\n    oof = []\n    \n    for train_index, test_index in folds.split(y):\n        \n        X_train, X_val = X.iloc[train_index], X.iloc[test_index]\n        y_train, y_val = y[train_index], y[test_index]\n        \n        #make and fit simple classifier model\n        model = lgb.LGBMRegressor(n_estimators = 50000, n_jobs = -1)\n        model.fit(X_train, y_train,\n          eval_set=[(X_train, y_train), (X_val, y_val)], eval_metric='auc',\n          verbose=0, early_stopping_rounds=500)\n        \n        models.append(model)\n        oof.append(roc_auc_score(y_val, model.predict(X_val)))\n        feature_importance += model.feature_importances_/N_FOLD\n    \n    preds = np.zeros(len(X))\n    test_preds = np.zeros(len(X_test))\n    \n    #predictions for train and test with these models\n    #column names lie halfway through the range, i.e. for y==1 at 3s to 4s, column name is '3.5'\n    for model in models:\n        preds += model.predict(X)\n        preds /= len(models)\n        df[str(threshold + sep/2)] = preds\n        \n        test_preds += model.predict(X_test)\n        test_preds /= len(models)\n        test_df[str(threshold + sep/2)] = test_preds\n    \n    #return the AUC for the combined classifiers at the target range            \n    return(np.asarray(oof).mean())","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"I'll be trying with a range of 0.5 seconds but this is easily modifiable for your own work. "},{"metadata":{"trusted":true},"cell_type":"code","source":"auc_af = []\nthresh = []\n#round maximum value down to smallest multiple of 0.5\nMAX_Y = np.floor(y_train.max()*2)/2\n\nfor i in np.arange(0, MAX_Y, SEP):\n    auc_af.append(ttf_classifier(i, X=train, X_test=test, df=train_preds, test_df=test_preds))\n    thresh.append(i)   ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"First let's take a look at the feature importances to see if anything is radically different for the classifiers:"},{"metadata":{"trusted":true},"cell_type":"code","source":"feature_importance /= len(auc_af)\nfeat_df = pd.DataFrame({'feature' : train.columns,\n                       'importance' : feature_importance}).sort_values('importance', ascending=False)\n\nplt.figure(figsize=(16, 32));\nsns.barplot(x='importance', y='feature', data=feat_df.iloc[:, :])\nplt.title('Mean Feature Importance')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"These look similar to the feature importances Andrew's Features data usually produce. \n\nHow well do the classifiers work at each 1s interval of `time_to_failure`? "},{"metadata":{"trusted":true},"cell_type":"code","source":"auc_af_df = pd.DataFrame({'threshold' : np.asarray(thresh) + SEP/2,\n                      'auc' : np.asarray(auc_af)})\n\nauc_af_df.plot(kind='line', x = 'threshold', y='auc', figsize=(15, 6))\nplt.title('AUC for predicting TTF lies in bin of width 1s')\nplt.xlabel('TTF')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"It looks like high TTF values are the easiest to identify. This dips around the central values, and the AUC climbs back up as TTF approaches 0, although it trails off at the end. At the extreme end of TTF there just aren't enough samples to reliably make predictions.\n\nMy current understanding of the data leads me to conclude:\n\n* Extremely low TTF, less than 1s, is harder to predict due to the gap between the slip/earthquake event and the measure of failure occuring. The acoustic profile of the system after the quake at this stage is similar to that at a high-TTF.\n* The 6-8s range is the most challenging for models due to the earthquake segments that contain 'micro-quakes' at this time\n\nNow we can evaluate how well the classifiers actually performed:"},{"metadata":{"trusted":true},"cell_type":"code","source":"#return column name with maximum probability for each row for train and tets\ntrain_preds['Pred'] = train_preds.apply(lambda x: x.argmax(), axis=1)\ntrain_preds['Pred'] = train_preds['Pred'].astype('float16')\n\ntest_preds['Pred'] = test_preds.apply(lambda x: x.argmax(), axis=1)\ntest_preds['Pred'] = test_preds['Pred'].astype('float16')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_preds['idx'] = train_preds.index\ntrain_preds['ttf'] = y_train\nax = train_preds.plot(kind='scatter',x='idx', y='Pred', figsize=(15, 6), s=1.5, color='b')\ntrain_preds.plot(kind='scatter',x='idx',y='ttf', figsize=(15, 6), s=0.75, color='r', ax=ax)\nplt.ylabel('TTF')\nplt.title('Predictions/Actual')\nplt.show()\n\nprint('MAE for classifier: ', mean_absolute_error(y_train, train_preds['Pred'].values))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The model evaluation graph looks awful! But the model itself is surprisingly (suspiciously) accurate. The graph above actually misrepresents how well the classifiers performed, since the accurate predictions are clumped together and the outliers are visually more obvious.\n\nMost competitors have found that their LB score was significantly better than their train-CV score. Sadly this is not the case for the classifiers, and you can expect an LB score of around 2 with this approach. Still, this method may have its uses. For example, it was very good at identifying the upper range of `time_to_failure`, and could be used in conjunction with your other models, which in my experience are more conservative at predicting either a very high or very low TTF. "},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"sub = pd.read_csv('../input/LANL-Earthquake-Prediction/sample_submission.csv')\nsub['time_to_failure'] = test_preds['Pred'].values\nsub.to_csv('utterly_dreadful_classifier_submission.csv', index=False)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### TTF Inequality Classifiers\n---\nWhat if instead we wanted to separate the data into two distinct ranges and run different models/algorithms on them? How accurately can we predict whether a segment is above or below a certain threshold? An even simpler variant of the above code can be used:"},{"metadata":{"trusted":true},"cell_type":"code","source":"def ttf_ineq_classifier(threshold, X=train, y=y_train):\n    \n    y = y >= threshold\n    models = []\n    oof = []\n    \n    for train_index, test_index in folds.split(y):\n        \n        X_train, X_val = X.iloc[train_index], X.iloc[test_index]\n        y_train, y_val = y[train_index], y[test_index]\n        \n        model = lgb.LGBMRegressor(n_estimators = 50000, n_jobs = -1)\n        model.fit(X_train, y_train,\n          eval_set=[(X_train, y_train), (X_val, y_val)], eval_metric='auc',\n          verbose=0, early_stopping_rounds=500)\n        models.append(model)\n        \n        oof.append(roc_auc_score(y_val, model.predict(X_val)))\n                           \n    return(np.asarray(oof).mean())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"SEP=0.1\nauc_ineq = []\nthresh = []\nfor i in np.arange(SEP, MAX_Y-1, SEP):\n    auc_ineq.append(ttf_ineq_classifier(i))\n    thresh.append(i)\n    \nauc_ineq_df = pd.DataFrame({'threshold' : np.asarray(thresh),\n                      'auc' : np.asarray(auc_ineq)})\n\nauc_ineq_df.plot(kind='line', x = 'threshold', y='auc', figsize=(15, 6))\nplt.title('AUC for predicting TTF > threshold')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"This corroborates the earlier findings that TTF is easiest to classify at its highest values. This also provides a more nuanced view of the TTF < 2s zone. While the basic classifiers found it hardest to identify TTF smaller than 2, once *extremely* close to failure, it is in fact a lot easier to identify. Strangely, the wrapper for LGBMClassifier worked far worse than the wrapper for LGBMRegressor. With the objective set to 'auc' this distinction largely vanishes, so this may be a product of the different hyperparameters they start with. \n\nHopefully these results will be of some use to you - remember that this kernel only used Andrew's Features so your own engineered data may yield more helpful results. Good luck!"}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.4","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}