{"cells":[{"metadata":{"_uuid":"fbec59384e3db10ebfb9ec0a1aa5a90b443e04de"},"cell_type":"markdown","source":"**Flavors of Physics**\n\nThis is a problem that is not easy to understand for someone who is not into Quantum Physics. Though a bit intriguing, the problem seems to be very interesting. \n\nIn simple terms, as I could undertand, the scientists are trying to observe the decay of tau to three muons - a phenomenon that is very rare.  The null hypothesis is that this decay does not happen and the idea here is to see that the decay happens more than what is currently known.\n\nThe dataset consists of real and simulated events. Real events (or Background events or those with Signal value = 0) and Simulated Events (or Signal Events or those with Signal value = 1) are part of the dataset. The analysis is to be done on the entire dataset to classify the data based on various attributes. \n\nBefore the algorithm is evaluated (based on AUC metric), the algorithm needs to pass through Agreement Test and Correlation Test. The code for these tests is provided along with the dataset."},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load in \n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# Input data files are available in the \"../input/\" directory.\n# For example, running this (by clicking run or pressing Shift+Enter) will list the files in the input directory\n\nimport os\nprint(os.listdir(\"../input/\"))\n\n# Any results you write to the current directory are saved as output.","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":true},"cell_type":"code","source":"import pandas as pd\nimport numpy\nimport seaborn as sns\nfrom sklearn.metrics import roc_curve, auc\nfrom sklearn.ensemble import GradientBoostingClassifier\nimport featuretools as ft\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"bc51d743f91c18093183c59fad9b60148a1a08ee"},"cell_type":"code","source":"def __rolling_window(data, window_size):\n    \"\"\"\n    Rolling window: take window with definite size through the array\n\n    :param data: array-like\n    :param window_size: size\n    :return: the sequence of windows\n\n    Example: data = array(1, 2, 3, 4, 5, 6), window_size = 4\n        Then this function return array(array(1, 2, 3, 4), array(2, 3, 4, 5), array(3, 4, 5, 6))\n    \"\"\"\n    shape = data.shape[:-1] + (data.shape[-1] - window_size + 1, window_size)\n    strides = data.strides + (data.strides[-1],)\n    return numpy.lib.stride_tricks.as_strided(data, shape=shape, strides=strides)\n\n\ndef __cvm(subindices, total_events):\n    \"\"\"\n    Compute Cramer-von Mises metric.\n    Compared two distributions, where first is subset of second one.\n    Assuming that second is ordered by ascending\n\n    :param subindices: indices of events which will be associated with the first distribution\n    :param total_events: count of events in the second distribution\n    :return: cvm metric\n    \"\"\"\n    target_distribution = numpy.arange(1, total_events + 1, dtype='float') / total_events\n    subarray_distribution = numpy.cumsum(numpy.bincount(subindices, minlength=total_events), dtype='float')\n    subarray_distribution /= 1.0 * subarray_distribution[-1]\n    return numpy.mean((target_distribution - subarray_distribution) ** 2)\n\n\ndef compute_cvm(predictions, masses, n_neighbours=200, step=50):\n    \"\"\"\n    Computing Cramer-von Mises (cvm) metric on background events: take average of cvms calculated for each mass bin.\n    In each mass bin global prediction's cdf is compared to prediction's cdf in mass bin.\n\n    :param predictions: array-like, predictions\n    :param masses: array-like, in case of Kaggle tau23mu this is reconstructed mass\n    :param n_neighbours: count of neighbours for event to define mass bin\n    :param step: step through sorted mass-array to define next center of bin\n    :return: average cvm value\n    \"\"\"\n    predictions = numpy.array(predictions)\n    masses = numpy.array(masses)\n    assert len(predictions) == len(masses)\n\n    # First, reorder by masses\n    predictions = predictions[numpy.argsort(masses)]\n\n    # Second, replace probabilities with order of probability among other events\n    predictions = numpy.argsort(numpy.argsort(predictions, kind='mergesort'), kind='mergesort')\n\n    # Now, each window forms a group, and we can compute contribution of each group to CvM\n    cvms = []\n    for window in __rolling_window(predictions, window_size=n_neighbours)[::step]:\n        cvms.append(__cvm(subindices=window, total_events=len(predictions)))\n    return numpy.mean(cvms)\n\n\ndef __roc_curve_splitted(data_zero, data_one, sample_weights_zero, sample_weights_one):\n    \"\"\"\n    Compute roc curve\n\n    :param data_zero: 0-labeled data\n    :param data_one:  1-labeled data\n    :param sample_weights_zero: weights for 0-labeled data\n    :param sample_weights_one:  weights for 1-labeled data\n    :return: roc curve\n    \"\"\"\n    labels = [0] * len(data_zero) + [1] * len(data_one)\n    weights = numpy.concatenate([sample_weights_zero, sample_weights_one])\n    data_all = numpy.concatenate([data_zero, data_one])\n    fpr, tpr, _ = roc_curve(labels, data_all, sample_weight=weights)\n    return fpr, tpr\n\n\ndef compute_ks(data_prediction, mc_prediction, weights_data, weights_mc):\n    \"\"\"\n    Compute Kolmogorov-Smirnov (ks) distance between real data predictions cdf and Monte Carlo one.\n\n    :param data_prediction: array-like, real data predictions\n    :param mc_prediction: array-like, Monte Carlo data predictions\n    :param weights_data: array-like, real data weights\n    :param weights_mc: array-like, Monte Carlo weights\n    :return: ks value\n    \"\"\"\n    assert len(data_prediction) == len(weights_data), 'Data length and weight one must be the same'\n    assert len(mc_prediction) == len(weights_mc), 'Data length and weight one must be the same'\n\n    data_prediction, mc_prediction = numpy.array(data_prediction), numpy.array(mc_prediction)\n    weights_data, weights_mc = numpy.array(weights_data), numpy.array(weights_mc)\n\n    assert numpy.all(data_prediction >= 0.) and numpy.all(data_prediction <= 1.), 'Data predictions are out of range [0, 1]'\n    assert numpy.all(mc_prediction >= 0.) and numpy.all(mc_prediction <= 1.), 'MC predictions are out of range [0, 1]'\n\n    weights_data /= numpy.sum(weights_data)\n    weights_mc /= numpy.sum(weights_mc)\n\n    fpr, tpr = __roc_curve_splitted(data_prediction, mc_prediction, weights_data, weights_mc)\n\n    Dnm = numpy.max(numpy.abs(fpr - tpr))\n    return Dnm\n\n\ndef roc_auc_truncated(labels, predictions, tpr_thresholds=(0.2, 0.4, 0.6, 0.8),\n                      roc_weights=(4, 3, 2, 1, 0)):\n    \"\"\"\n    Compute weighted area under ROC curve.\n\n    :param labels: array-like, true labels\n    :param predictions: array-like, predictions\n    :param tpr_thresholds: array-like, true positive rate thresholds delimiting the ROC segments\n    :param roc_weights: array-like, weights for true positive rate segments\n    :return: weighted AUC\n    \"\"\"\n    assert numpy.all(predictions >= 0.) and numpy.all(predictions <= 1.), 'Data predictions are out of range [0, 1]'\n    assert len(tpr_thresholds) + 1 == len(roc_weights), 'Incompatible lengths of thresholds and weights'\n    fpr, tpr, _ = roc_curve(labels, predictions)\n    area = 0.\n    tpr_thresholds = [0.] + list(tpr_thresholds) + [1.]\n    for index in range(1, len(tpr_thresholds)):\n        tpr_cut = numpy.minimum(tpr, tpr_thresholds[index])\n        tpr_previous = numpy.minimum(tpr, tpr_thresholds[index - 1])\n        area += roc_weights[index - 1] * (auc(fpr, tpr_cut, reorder=True) - auc(fpr, tpr_previous, reorder=True))\n    tpr_thresholds = numpy.array(tpr_thresholds)\n    # roc auc normalization to be 1 for an ideal classifier\n    area /= numpy.sum((tpr_thresholds[1:] - tpr_thresholds[:-1]) * numpy.array(roc_weights))\n    return area\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7f1cd5a6e6df4942ca04eb3b1d490c8d288d02a4"},"cell_type":"code","source":"folder = \"../input/\"\ntrain = pd.read_csv(folder+'training.csv', index_col='id')\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5cdaacb0dec61c4723199a1edf3e85962fb21a9c"},"cell_type":"code","source":"train.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a87557c500edb5421b97585e1531824aa1e427dc"},"cell_type":"code","source":"test = pd.read_csv(folder + 'test.csv', index_col='id')\ntest.shape","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"4f1c4b76d138a67d47cd0ab9700f7c1275c03d23"},"cell_type":"markdown","source":"Let's try to visualise the data - it always helps to see if there is some apparent relationship or anomaly in the dataset."},{"metadata":{"trusted":true,"_uuid":"55af4721b23b74e64f40e30b66832385e7b40130"},"cell_type":"code","source":"# visualize the relationship between the features and the response using scatterplots\nsns.pairplot(train, x_vars=['LifeTime','dira','FlightDistance','FlightDistanceError'], y_vars='signal', size=7, aspect=0.7)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4636d3ca78f412d43437ca0fa7801607cadad5cd"},"cell_type":"code","source":"sns.pairplot(train, x_vars=['LifeTime','dira','FlightDistance'], y_vars='FlightDistanceError', size=7, aspect=0.7)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"67efffa7d614a224e58ff9d478b191009a5527d5"},"cell_type":"markdown","source":"Since the test dataset does not have 4 features, it is prudent to remove these from the training dataset for the analysis."},{"metadata":{"trusted":true,"_uuid":"8fdc9605fde85ea2b9ae506bf5667ae0600cdba8","scrolled":true},"cell_type":"code","source":"removeFeatures = ('signal','production','mass','min_ANNmuon')\ntrTarget = train.loc[:,'signal']\ntrTarget.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"29cb876b9327f2a626c7d67ac144b7f79a70386e"},"cell_type":"code","source":"trainNew = train.drop(columns=['signal','production','mass','min_ANNmuon'])\ntrainNew.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"9741bc03490fdcce2069ffcc073ecd997c394293"},"cell_type":"code","source":"trainNew.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"dbf11fe4ad268313b6546d6758629f7cd87c8985"},"cell_type":"code","source":"trainNew.tail(1)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5e0cc33eb0c2e95534c108e9318387611dfa37e5"},"cell_type":"code","source":"test.head(1)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"2f747258ab62a41bcacf86fe007001563964a918"},"cell_type":"markdown","source":"Let's try to add some features to the dataset. \n\nWe can add features by introducing some categorical variable (if binary, we do not need OneHotEncoding). \nWe can also use FeatureTools library to automatically generate features."},{"metadata":{"trusted":true,"_uuid":"8bbedf3e20bfd0eaa3710bd7bb7d95e99dbcbf1c"},"cell_type":"code","source":"sns.boxplot( y=trainNew[\"pt\"] )\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7a4e5a77087011f8d42f1725825afd6cece40be7"},"cell_type":"code","source":"trainNew['pt'].describe()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"183c2294b162e0d481a5882f9433ab906f8d020d"},"cell_type":"code","source":"combi = pd.concat([trainNew, test])\ncombi.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7247fc2542220031b6a4ce712d66ade5504bde33"},"cell_type":"code","source":"def convertNum2Bin(x):\n    if 0 < x <= 5028:\n        return 0\n    return 1\n\ncombi['ptBin'] = combi['pt'].apply(convertNum2Bin)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1bbd48ae8244dc78398a3e1117f292ecbd44ec51"},"cell_type":"code","source":"trainNew['LifeTime'].describe()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"0ce512f6a8ad99ce1c519b5ca365f0f11d230ab1"},"cell_type":"code","source":"def convertLT2Bin(x):\n    if 0 < x <= 0.001255:\n        return 0\n    return 1\n\ncombi['LifeTimeBin'] = combi['LifeTime'].apply(convertLT2Bin)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7819782cc8bab188c2941d9350a7d0784f795d99"},"cell_type":"code","source":"trainNew['FlightDistance'].describe()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"3149af1eb96a1d2aa9c3601490679539f3718523"},"cell_type":"code","source":"def convertFD2Bin(x):\n    if 0 < x <= 15.154:\n        return 0\n    return 1\n\ncombi['FDBin'] = combi['FlightDistance'].apply(convertFD2Bin)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a9d00f0c6e66b20e9e0f95afa8278d9261b8cf42"},"cell_type":"code","source":"trainNew['FlightDistanceError'].describe()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7495d230facfc619654bbb46e66775402e182a92"},"cell_type":"code","source":"def convertFDE2Bin(x):\n    if 0 < x <= 0.5018:\n        return 0\n    return 1\n\ncombi['FDEBin'] = combi['FlightDistanceError'].apply(convertFDE2Bin)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"bdc7f99ca414d544fd3f68fab8a5368db5623cf6"},"cell_type":"code","source":"combi.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1ef1f2d0ecde962d864a03f7d572ae430d83c1fa"},"cell_type":"code","source":"\n# Make an entityset and add the entity\n#es = ft.EntitySet(id = 'particleData')\n#es.entity_from_dataframe(entity_id = 'TrainAndTest', dataframe = combi, \n#                         make_index = True, index = 'index')\n\n# Run deep feature synthesis with transformation primitives\n#feature_matrix, feature_defs = ft.dfs(entityset = es, target_entity = 'TrainAndTest', max_depth = 2, verbose = 1,\n                                      agg_primitives = ['count','mean'])\n#                                      trans_primitives = ['add', 'multiply'])\n\n#feature_matrix.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"da225096398d3d623b87682341020566cd9d6dd7"},"cell_type":"code","source":"#feature_matrix.columns\n#feature_matrix.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a243a7561c8e614f4bd2a6970515936b3667cb59"},"cell_type":"code","source":"#combi.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d98c1c5fa37064935f63efbceb6c4c8fc161389c"},"cell_type":"code","source":"#combi['id']=combi.reset_index().index\n#combi.iloc[67552:67554,:]\ntrainNew = combi.iloc[:67553,:]    # or should we use feature_matrix instead of combi?\ntest= combi.iloc[67553:,:]\ntrainNew.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4e63be990aba8beee8349a3515382e148c0b566a"},"cell_type":"code","source":"test.shape","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"26a12a8d0bfcf57dc9f8c997d06fde5e45e075d4"},"cell_type":"markdown","source":"Before we proceed, let's Scale the data so that all features contribute equally to the prediction. We do not want some features to have more weight compared to others just becasue of a difference in the units of measurement."},{"metadata":{"trusted":true,"_uuid":"7cd57007839fe879bcac6a726a3316178b1c8f9b"},"cell_type":"code","source":"from sklearn.preprocessing import StandardScaler\n\nscaler = StandardScaler()\ntrainScaled = scaler.fit_transform(trainNew)\ntrainScaled[0:2,]\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7a0d7ea0a1770b50f22831c9acb2369a5cb0af5a"},"cell_type":"code","source":"testScaled = scaler.transform(test)\ntestScaled[0:2,]\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"30013af37ec31ce8979ab9e924e36745de6f53ef"},"cell_type":"markdown","source":"It might be a good idea to do a Principal Component Analysis of the data to avoid features that are highly correlated. Also, we will go for features that represent 95% of the variance to avoid overfitting of the data."},{"metadata":{"trusted":true,"_uuid":"063711ce73c899835d2308f32e9dda396845c4f1"},"cell_type":"code","source":"from sklearn.decomposition import PCA\n\npca = PCA().fit(trainScaled)\ntrain_pca = pca.transform(trainScaled)\ntest_pca = pca.transform(testScaled)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"9d64db3fab9c7539900b95a08c459ddb6b899519"},"cell_type":"code","source":"def pca_summary(pca, standardized_data, out=True):\n    names = [\"PC\"+str(i) for i in range(1, len(pca.explained_variance_ratio_)+1)]\n    a = list(np.std(pca.transform(standardized_data), axis=0))\n    b = list(pca.explained_variance_ratio_)\n    c = [np.sum(pca.explained_variance_ratio_[:i]) for i in range(1, len(pca.explained_variance_ratio_)+1)]\n    columns = pd.MultiIndex.from_tuples([(\"sdev\", \"Standard deviation\"), (\"varprop\", \"Proportion of Variance\"), (\"cumprop\", \"Cumulative Proportion\")])\n    summary = pd.DataFrame(list(zip(a, b, c)), index=names, columns=columns)\n    \n    if out:\n        print(\"Importance of components:\")\n        display(summary)\n    return summary","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"649742d26bcaf1ecb0fa5d162e710fd6dd13ba12"},"cell_type":"code","source":"summary = pca_summary(pca, train_pca)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e1f9d94e67327e12987c18127c2d8945e581c54f"},"cell_type":"code","source":"import matplotlib.pyplot as plt\n\ndef screeplot(pca, standardized_values):\n    y = np.std(pca.transform(standardized_values), axis=0)**2\n    x = np.arange(len(y)) + 1\n    plt.plot(x, y, \"o-\")\n    plt.xticks(x, [\"Comp.\"+str(i) for i in x], rotation=60)\n    plt.ylabel(\"Variance\")\n    plt.show()\n\nscreeplot(pca, trainScaled)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"6a1151c66093f72d57f7a1a170883028217cf823"},"cell_type":"markdown","source":"First 31 PCA attributes give more than 95% variance. Let's go ahead with 31 PCA Attributes for our analysis."},{"metadata":{"trusted":true,"scrolled":true,"_uuid":"72bab283dc1ac9aae64d8928b1449b8c99845241"},"cell_type":"code","source":"pcatrain = pd.DataFrame(train_pca[:,0:31])\npcatrain.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"3252d2572c05bc506d4dbd0e9c558979d6dcfe94"},"cell_type":"code","source":"pcatest = pd.DataFrame(test_pca[:,0:31])\npcatest.head()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"dae8f61cc20f8df1e3ec468b6b6b06bbc7bb18b5"},"cell_type":"markdown","source":"We will try various models - Gradient Boosting Classifier, Random Forest Classifier and Extreme Gradient Boosting (XGB) Classifier. \n\nWe can try Parameter Tuning using Grid Cross Validation to fine tune the model. "},{"metadata":{"_uuid":"c97cc59306dcfaca1e39aa50804fa33a0c4a94b2"},"cell_type":"markdown","source":"**Baseline Training (Gradient Boosting Classifier)**"},{"metadata":{"trusted":true,"_uuid":"4d0533d8b8f21208fa14796010a0f4731472b49f"},"cell_type":"code","source":"baseline = GradientBoostingClassifier(n_estimators=40, learning_rate=0.01, subsample=0.7,\n                                      min_samples_leaf=10, max_depth=7, random_state=11)\nbaseline.fit(pcatrain, trTarget)\n#baseline.fit(train[variables], train['signal'])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"7160c8b706a4af3ad47a04f020cb0725dfc42f14"},"cell_type":"markdown","source":"**Random Forest Classifier (with Parameter Tuning)**"},{"metadata":{"trusted":true,"_uuid":"404186808a772ac353efdd8694d4c24d84664105"},"cell_type":"code","source":"from sklearn.ensemble import RandomForestClassifier\nfrom sklearn.model_selection import GridSearchCV\n#pcatrain, pcatest, trTarget - already available\n\n# Create the parameter grid \nparam_grid = {\n    'bootstrap': [False],\n#    'max_depth': [100, 200],       ------ The RF Classifier failed the Agreement Test based on these parameters. So, leaving this step out for faster re-run of code.\n    'max_features': ['sqrt'],\n#    'min_samples_leaf': [1, 3],\n#    'n_estimators': [200, 400]\n}\n# Create a RF Classifier model\nrf = RandomForestClassifier(random_state = 42)\n# Instantiate the grid search model\n#grid_search = GridSearchCV(estimator = rf, param_grid = param_grid, cv = 2, n_jobs = -1, verbose = 2) --- CV works but takes time so doing it without CV\ngrid_search = GridSearchCV(estimator = rf, param_grid = param_grid, n_jobs = -1, verbose = 2)\n# Fit the grid search to the data\ngrid_search.fit(pcatrain, trTarget)\ngrid_search.best_params_\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a8995975b1d4b9d613420ae5df3063a44ab0033e"},"cell_type":"code","source":"best_gridRFCModel = grid_search.best_estimator_","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"3183de752fc59c394655bbd9cba498af774e1e08"},"cell_type":"markdown","source":" **XGB Classifier (First set of parameters)**"},{"metadata":{"trusted":true,"_uuid":"5b3090da8a4457e33b59188638df2a0df05a47fe"},"cell_type":"code","source":"import xgboost as xgb\n\n\nclfXGB = xgb.XGBClassifier(learning_rate =0.1, n_estimators=1000, max_depth=5,random_state = 42)\nclfXGB.fit(pcatrain, trTarget)\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"6c4d11b28cafa3a388944e720c1145784db929e9"},"cell_type":"markdown","source":"**XGB Classifier (Second set of parameters)**"},{"metadata":{"trusted":true,"_uuid":"7e1b8b23f5c81405afffa3cbfd4a766f47ed0916"},"cell_type":"code","source":"clfXGB2 = xgb.XGBClassifier(n_estimators=40, learning_rate=0.01, subsample=0.7,min_samples_leaf=10, max_depth=7, random_state=11)\nclfXGB2.fit(pcatrain, trTarget)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a4866266195739139e6ae7221c98294b7c9bd500"},"cell_type":"markdown","source":" **XGB Classifier (Third set of parameters)**"},{"metadata":{"trusted":true,"_uuid":"00d5075777a08275d8201cfa2e28d9c4649470c0"},"cell_type":"code","source":"clfXGB3 = xgb.XGBClassifier(n_estimators=50, max_depth = 9, learning_rate=0.01, subsample=0.75, random_state=11)\nclfXGB3.fit(pcatrain, trTarget)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"7d4d4fd91707abdb5a36df1cbb08a7c0efe99a3f"},"cell_type":"markdown","source":"# Check agreement test"},{"metadata":{"trusted":true,"_uuid":"3b400e987a7906cc203d931c6ec0600b8e146682"},"cell_type":"code","source":"check_agreement = pd.read_csv(folder + 'check_agreement.csv', index_col='id')\ncheck_agreement.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"398ec4a4faefc6cb1af19780d4d731324805133c"},"cell_type":"code","source":"check_agreement.head(2)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"9124e535b0dc509d00115590384ef524bab3f9fb"},"cell_type":"markdown","source":"Check Agreement Dataframe has weight and signal in addition to the 50 attributes given in Train dataset."},{"metadata":{"trusted":true,"_uuid":"1199083cbdd8e4340f57aacf678d8517880db329"},"cell_type":"code","source":"check_agreement_sig_wt = check_agreement.loc[:,['signal','weight']]\ncheck_agreement_sig_wt.head(2)\n\ncheck_agreement_var = check_agreement.drop(columns=['signal','weight']) \ncheck_agreement_var.head(2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"563088782da117399a49d86c4a3c9c9ed948fd3e"},"cell_type":"code","source":"check_agreement_var['ptBin'] = check_agreement_var['pt'].apply(convertNum2Bin)\ncheck_agreement_var['LifeTimeBin'] = check_agreement_var['LifeTime'].apply(convertLT2Bin)\ncheck_agreement_var['FDBin'] = check_agreement_var['FlightDistance'].apply(convertFD2Bin)\ncheck_agreement_var['FDEBin'] = check_agreement_var['FlightDistanceError'].apply(convertFDE2Bin)\n\nagreementScaled = scaler.transform(check_agreement_var)\nagreement_pca = pca.transform(agreementScaled)\npcaAgreement = pd.DataFrame(agreement_pca[:,0:31])\n\npcaAgreement.head(2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"f45970d0a4f6111e5d1d40971704507894f40b2b"},"cell_type":"code","source":"#agreement_probs = baseline.predict_proba(check_agreement[variables])[:, 1]\nagreement_probs = baseline.predict_proba(pcaAgreement)[:, 1]\n\nks = compute_ks(\n    agreement_probs[check_agreement_sig_wt['signal'].values == 0],\n    agreement_probs[check_agreement_sig_wt['signal'].values == 1],\n    check_agreement[check_agreement_sig_wt['signal'] == 0]['weight'].values,\n    check_agreement[check_agreement_sig_wt['signal'] == 1]['weight'].values)\nprint ('KS metric', ks, ks < 0.09)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"7073c4dca1c1f04a14cecf4adc245a8dad88c8c2"},"cell_type":"markdown","source":"Let's do the Agreement Test on Random Forest Classifier too."},{"metadata":{"trusted":true,"_uuid":"72e83ad1616b85558d74a1f6f667301e5597f4fe"},"cell_type":"code","source":"agreement_probs = best_gridRFCModel.predict_proba(pcaAgreement)[:, 1]\n\nks = compute_ks(\n    agreement_probs[check_agreement_sig_wt['signal'].values == 0],\n    agreement_probs[check_agreement_sig_wt['signal'].values == 1],\n    check_agreement[check_agreement_sig_wt['signal'] == 0]['weight'].values,\n    check_agreement[check_agreement_sig_wt['signal'] == 1]['weight'].values)\nprint ('KS metric', ks, ks < 0.09)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"50907da24370d94bcc9ceaab56df62134db56f03"},"cell_type":"markdown","source":"Random Forest Classifier test passed the test. Let's do it on XGB Classifier for 3 different sets of parameters."},{"metadata":{"trusted":true,"_uuid":"0eae127d063093bbcd8e0f0a4add51a837bfd2d6"},"cell_type":"code","source":"agreement_probs = clfXGB.predict_proba(pcaAgreement)[:, 1]\n\nks = compute_ks(\n    agreement_probs[check_agreement_sig_wt['signal'].values == 0],\n    agreement_probs[check_agreement_sig_wt['signal'].values == 1],\n    check_agreement[check_agreement_sig_wt['signal'] == 0]['weight'].values,\n    check_agreement[check_agreement_sig_wt['signal'] == 1]['weight'].values)\nprint ('KS metric', ks, ks < 0.09)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"90ef6624f9fd0796ee17587ea50c7c79ad507289"},"cell_type":"code","source":"agreement_probs = clfXGB2.predict_proba(pcaAgreement)[:, 1]\n\nks = compute_ks(\n    agreement_probs[check_agreement_sig_wt['signal'].values == 0],\n    agreement_probs[check_agreement_sig_wt['signal'].values == 1],\n    check_agreement[check_agreement_sig_wt['signal'] == 0]['weight'].values,\n    check_agreement[check_agreement_sig_wt['signal'] == 1]['weight'].values)\nprint ('KS metric', ks, ks < 0.09)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"8706bf05ffb933092650e2a24c91e1181f1cf666"},"cell_type":"code","source":"agreement_probs = clfXGB3.predict_proba(pcaAgreement)[:, 1]\n\nks = compute_ks(\n    agreement_probs[check_agreement_sig_wt['signal'].values == 0],\n    agreement_probs[check_agreement_sig_wt['signal'].values == 1],\n    check_agreement[check_agreement_sig_wt['signal'] == 0]['weight'].values,\n    check_agreement[check_agreement_sig_wt['signal'] == 1]['weight'].values)\nprint ('KS metric', ks, ks < 0.09)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"f3e7de2f73fdffec397043cd074f7504f613daee"},"cell_type":"markdown","source":"Out of the above 5 models, Baseline model (Gradient Boosting), Random Forest Classifier, XGB (with second set of parameters) and XGB (with third set of parameters) passed the Agreement Test. We will proceed with these 4 models for the next steps."},{"metadata":{"_uuid":"69663d276c5e4927b4240b81878251ff1ccf0cf5"},"cell_type":"markdown","source":"# Check correlation test\n"},{"metadata":{"trusted":true,"_uuid":"13b043b17e60ec4e079d9d33176c5a91220b62e7"},"cell_type":"code","source":"check_correlation = pd.read_csv(folder + 'check_correlation.csv', index_col='id')\ncheck_correlation.shape\ncheck_correlation.head(1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a02380f9ba5a57b4d4c7ffe4de7199a82800efef"},"cell_type":"code","source":"check_correlation_mass = check_correlation.loc[:,'mass']\ncheck_correlation_mass.head(2)\n\ncheck_correlation_var = check_correlation.drop('mass',axis = 1) \ncheck_correlation_var.head(2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"37f5767ece04a91eede2e4a9088c613814992366"},"cell_type":"code","source":"check_correlation_var.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"241ea283ba14f3b367703dbc32fc0eff7aaa0846"},"cell_type":"code","source":"check_correlation_var['ptBin'] = check_correlation_var['pt'].apply(convertNum2Bin)\ncheck_correlation_var['LifeTimeBin'] = check_correlation_var['LifeTime'].apply(convertLT2Bin)\ncheck_correlation_var['FDBin'] = check_correlation_var['FlightDistance'].apply(convertFD2Bin)\ncheck_correlation_var['FDEBin'] = check_correlation_var['FlightDistanceError'].apply(convertFDE2Bin)\n\ncorrelationScaled = scaler.transform(check_correlation_var)\ncorrelation_pca = pca.transform(correlationScaled)\npcaCorrelation = pd.DataFrame(correlation_pca[:,0:31])\n\npcaCorrelation.head(2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"99af0f3d803626ab837d1e155abfdf523a72987d"},"cell_type":"code","source":"correlation_probs = baseline.predict_proba(pcaCorrelation)[:, 1]\ncvm = compute_cvm(correlation_probs, check_correlation_mass)\nprint ('CvM metric', cvm, cvm < 0.002)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"101d1a2a24aef7582fb2315ea56b45f9ecff7d1a"},"cell_type":"code","source":"correlation_probs = best_gridRFCModel.predict_proba(pcaCorrelation)[:, 1]\ncvm = compute_cvm(correlation_probs, check_correlation_mass)\nprint ('CvM metric', cvm, cvm < 0.002)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"636d606a9c203781249ff1671d3edc088d18a5b9"},"cell_type":"code","source":"correlation_probs = clfXGB2.predict_proba(pcaCorrelation)[:, 1]\ncvm = compute_cvm(correlation_probs, check_correlation_mass)\nprint ('CvM metric', cvm, cvm < 0.002)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a97cb14dc928e37028b4b2a50a91c348a2fef6d5"},"cell_type":"code","source":"correlation_probs = clfXGB3.predict_proba(pcaCorrelation)[:, 1]\ncvm = compute_cvm(correlation_probs, check_correlation_mass)\nprint ('CvM metric', cvm, cvm < 0.002)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"27221caa6290d5a8902f1140816f0641d3e60c12"},"cell_type":"markdown","source":"So, only 3 of the models - passed the Correlation Test."},{"metadata":{"_uuid":"0ef1119873c771eeb6159632718230bc328e4a33"},"cell_type":"markdown","source":"# Compute weighted AUC on the training data with min_ANNmuon > 0.4\n"},{"metadata":{"trusted":true,"_uuid":"ac923b6ee1167903695532b76282bb554a379a3d"},"cell_type":"code","source":"train_eval = train[train['min_ANNmuon'] > 0.4]\ntrainEvalSignal = train_eval['signal']\ntrainEvalTruncated = train_eval.drop(columns=['signal','production','mass','min_ANNmuon'])\n\ntrainEvalTruncated['ptBin'] = trainEvalTruncated['pt'].apply(convertNum2Bin)\ntrainEvalTruncated['LifeTimeBin'] = trainEvalTruncated['LifeTime'].apply(convertLT2Bin)\ntrainEvalTruncated['FDBin'] = trainEvalTruncated['FlightDistance'].apply(convertFD2Bin)\ntrainEvalTruncated['FDEBin'] = trainEvalTruncated['FlightDistanceError'].apply(convertFDE2Bin)\n\ntrainEvalScaled = scaler.transform(trainEvalTruncated)\ntrainEval_pca = pca.transform(trainEvalScaled)\npcaTrainEval = pd.DataFrame(trainEval_pca[:,0:31])\n\npcaTrainEval.head(2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c83321c49c0e27efe97690c1a612956299ae2c91"},"cell_type":"code","source":"train_probs = baseline.predict_proba(pcaTrainEval)[:, 1]\nAUC = roc_auc_truncated(trainEvalSignal, train_probs)\nprint ('AUC', AUC)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7b89574ac6da858a7f7429d440745f7df8533575"},"cell_type":"code","source":"train_probs = clfXGB2.predict_proba(pcaTrainEval)[:, 1]\nAUC = roc_auc_truncated(trainEvalSignal, train_probs)\nprint ('AUC', AUC)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"b6c2062bfcdf07baedbefae7ed359780887dddf8"},"cell_type":"code","source":"train_probs = clfXGB3.predict_proba(pcaTrainEval)[:, 1]\nAUC = roc_auc_truncated(trainEvalSignal, train_probs)\nprint ('AUC', AUC)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a5a8796c794cff097b469786ab310f5f16133725"},"cell_type":"markdown","source":"AUC Metric is best for XGB with the third set of parameters. Let's create file for Kaggle."},{"metadata":{"_uuid":"85b103da9cbf8433ed8d042b27d809fbbc2a8ec3"},"cell_type":"markdown","source":"# Predict test, create file for kaggle"},{"metadata":{"trusted":true,"_uuid":"a1fc3fd359f50a403adad71deac0d8dca5191d4f"},"cell_type":"code","source":"result = pd.DataFrame({'id': test.index})\nresult['prediction'] = clfXGB3.predict_proba(pcatest)[:, 1]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"cf01f3915f6ddae44a363c9e1b6f3e0fe32ece8f"},"cell_type":"code","source":"sub = result.to_csv('Achal_baseline_5.csv', index=False, sep=',')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d8dee59b4e04b58fc6a5a6096423d86011389840"},"cell_type":"markdown","source":""}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.6","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}