{"cells":[{"metadata":{"_uuid":"8c3751874ef1429d607e9f23dee0e0024b166c65"},"cell_type":"markdown","source":"# Intro\nThis is my first kaggle competition and I have really enjoyed it so far.\nThe following kernel has been inspired by:\n+ https://www.kaggle.com/ogrellier/plasticc-in-a-kernel-meta-and-data\n+ https://www.kaggle.com/c/PLAsTiCC-2018/discussion/70908\n+ https://www.kaggle.com/meaninglesslives/simple-neural-net-for-time-series-classification\n\nStill haven't figured out the feature hidden in here https://www.kaggle.com/c/PLAsTiCC-2018/discussion/70725#416740\n\nA big thanks to the kaggle community that makes this competition so enjoyable."},{"metadata":{"_uuid":"609006646970a83c89ff8a332459ddceb2e02e4c"},"cell_type":"markdown","source":"# Loading Libraries\n\n\n"},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nfrom sklearn.metrics import log_loss\nfrom sklearn.model_selection import StratifiedKFold\nimport gc\nimport os\nimport matplotlib.pyplot as plt\nimport seaborn as sns \nimport lightgbm as lgb\nfrom catboost import Pool, CatBoostClassifier\nimport itertools\nimport pickle, gzip\nimport glob\nfrom sklearn.preprocessing import StandardScaler\nfrom tsfresh.feature_extraction import extract_features\nnp.warnings.filterwarnings('ignore')\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"7dbccfa992f29644a47a31341b2c68c6b42b835d"},"cell_type":"markdown","source":"# Extracting Features from train set\n+ Features from olivier's kernel \n+ very smart and powerful feature that is generously given here https://www.kaggle.com/c/PLAsTiCC-2018/discussion/69696#410538\n+ per passband features with tsfresh library. fft features added to capture periodicity https://www.kaggle.com/c/PLAsTiCC-2018/discussion/70346#415506"},{"metadata":{"trusted":true,"_uuid":"bc7cd924880e0c44cdb7130ca73ace3e0026de4c"},"cell_type":"code","source":"gc.enable()\n\ntrain = pd.read_csv('../input/training_set.csv')\n# Features to compute with tsfresh library. Fft coefficient is meant to capture periodicity\nfcp = {'fft_coefficient': [{'coeff': 0, 'attr': 'abs'},{'coeff': 1, 'attr': 'abs'}],'kurtosis' : None, 'skewness' : None}\n\ndef featurize(df):\n    df['flux_ratio_sq'] = np.power(df['flux'] / df['flux_err'], 2.0)\n    df['flux_by_flux_ratio_sq'] = df['flux'] * df['flux_ratio_sq']\n    # train[detected==1, mjd_diff:=max(mjd)-min(mjd), by=object_id]\n\n\n    aggs = {\n        'flux': ['min', 'max', 'mean', 'median', 'std','skew'],\n        'flux_err': ['min', 'max', 'mean', 'median', 'std','skew'],\n        'detected': ['mean'],\n        'flux_ratio_sq':['sum','skew'],\n        'flux_by_flux_ratio_sq':['sum','skew'],\n    }\n\n    agg_df = df.groupby('object_id').agg(aggs)\n    new_columns = [\n        k + '_' + agg for k in aggs.keys() for agg in aggs[k]\n    ]\n    agg_df.columns = new_columns\n    agg_df['flux_diff'] = agg_df['flux_max'] - agg_df['flux_min']\n    agg_df['flux_dif2'] = (agg_df['flux_max'] - agg_df['flux_min']) / agg_df['flux_mean']\n    agg_df['flux_w_mean'] = agg_df['flux_by_flux_ratio_sq_sum'] / agg_df['flux_ratio_sq_sum']\n    agg_df['flux_dif3'] = (agg_df['flux_max'] - agg_df['flux_min']) / agg_df['flux_w_mean']\n    # Add more features with \n    agg_df_ts = extract_features(df, column_id='object_id', column_sort='mjd', column_kind='passband', column_value = 'flux', default_fc_parameters = fcp, n_jobs=4)\n    # Add smart feature that is suggested here https://www.kaggle.com/c/PLAsTiCC-2018/discussion/69696#410538\n    # dt[detected==1, mjd_diff:=max(mjd)-min(mjd), by=object_id]\n    df_det = df[df['detected']==1].copy()\n\n    agg_df_mjd = extract_features(df_det, column_id='object_id', column_value = 'mjd', default_fc_parameters = {'maximum':None, 'minimum':None}, n_jobs=4)\n    agg_df_mjd['mjd_diff_det'] = agg_df_mjd['mjd__maximum'] - agg_df_mjd['mjd__minimum']\n    del agg_df_mjd['mjd__maximum'], agg_df_mjd['mjd__minimum']\n    agg_df_ts = pd.merge(agg_df_ts, agg_df_mjd, on = 'id')\n    # tsfresh returns a dataframe with an index name='id'\n    agg_df_ts.index.rename('object_id',inplace=True)\n    agg_df = pd.merge(agg_df, agg_df_ts, on='object_id')\n    return agg_df\n\nagg_train = featurize(train)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"7ed21e25fb5678c78bb19a8297bf83db00ccae01"},"cell_type":"markdown","source":"# Merging extracted features with meta data"},{"metadata":{"trusted":true,"_uuid":"696ee870837c23375bc7f0388de1b07510351123"},"cell_type":"code","source":"meta_train = pd.read_csv('../input/training_set_metadata.csv')\nmeta_train.head()\n\nfull_train = agg_train.reset_index().merge(\n    right=meta_train,\n    how='outer',\n    on='object_id'\n)\n\nif 'target' in full_train:\n    y = full_train['target']\n    del full_train['target']\nclasses = sorted(y.unique())\n\n# Taken from Giba's topic : https://www.kaggle.com/titericz\n# https://www.kaggle.com/c/PLAsTiCC-2018/discussion/67194\n# with Kyle Boone's post https://www.kaggle.com/kyleboone\nclass_weight = {\n    c: 1 for c in classes\n}\nfor c in [64, 15]:\n    class_weight[c] = 2\n\nprint('Unique classes : ', classes)\n\nif 'object_id' in full_train:\n    oof_df = full_train[['object_id']]\n    del full_train['object_id'], full_train['distmod'], full_train['hostgal_specz']\n    del full_train['ra'], full_train['decl'], full_train['gal_l'],full_train['gal_b'],full_train['ddf']\n    \ntrain_mean = full_train.mean(axis=0)\nfull_train.fillna(0, inplace=True)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"7aaadd7fb68f26be1aa78e0721332689365f2322"},"cell_type":"markdown","source":"# Class weights and loss functions"},{"metadata":{"trusted":true,"_uuid":"4c35b11a75cdc6aa926109ae876a1f2bb9f4078c"},"cell_type":"code","source":"# Compute weights\nw = y.value_counts()\nweights = {i : np.sum(w) / w[i] for i in w.index}\n\ndef multi_weighted_logloss(y_true, y_preds):\n    \"\"\"\n    @author olivier https://www.kaggle.com/ogrellier\n    multi logloss for PLAsTiCC challenge\n    \"\"\"\n    # class_weights taken from Giba's topic : https://www.kaggle.com/titericz\n    # https://www.kaggle.com/c/PLAsTiCC-2018/discussion/67194\n    # with Kyle Boone's post https://www.kaggle.com/kyleboone\n    classes = [6, 15, 16, 42, 52, 53, 62, 64, 65, 67, 88, 90, 92, 95]\n    class_weight = {6: 1, 15: 2, 16: 1, 42: 1, 52: 1, 53: 1, 62: 1, 64: 2, 65: 1, 67: 1, 88: 1, 90: 1, 92: 1, 95: 1}\n    if len(np.unique(y_true)) > 14:\n        classes.append(99)\n        class_weight[99] = 2\n    y_p = y_preds\n    # Trasform y_true in dummies\n    y_ohe = pd.get_dummies(y_true)\n    # Normalize rows and limit y_preds to 1e-15, 1-1e-15\n    y_p = np.clip(a=y_p, a_min=1e-15, a_max=1 - 1e-15)\n    # Transform to log\n    y_p_log = np.log(y_p)\n    # Get the log for ones, .values is used to drop the index of DataFrames\n    # Exclude class 99 for now, since there is no class99 in the training set\n    # we gave a special process for that class\n    y_log_ones = np.sum(y_ohe.values * y_p_log, axis=0)\n    # Get the number of positives for each class\n    nb_pos = y_ohe.sum(axis=0).values.astype(float)\n    # Weight average and divide by the number of positives\n    class_arr = np.array([class_weight[k] for k in sorted(class_weight.keys())])\n    y_w = y_log_ones * class_arr / nb_pos\n\n    loss = - np.sum(y_w) / np.sum(class_arr)\n    return loss\n\n\ndef lgb_multi_weighted_logloss(y_true, y_preds):\n    \"\"\"\n    @author olivier https://www.kaggle.com/ogrellier\n    multi logloss for PLAsTiCC challenge\n    \"\"\"\n    # class_weights taken from Giba's topic : https://www.kaggle.com/titericz\n    # https://www.kaggle.com/c/PLAsTiCC-2018/discussion/67194\n    # with Kyle Boone's post https://www.kaggle.com/kyleboone\n    classes = [6, 15, 16, 42, 52, 53, 62, 64, 65, 67, 88, 90, 92, 95]\n    class_weight = {6: 1, 15: 2, 16: 1, 42: 1, 52: 1, 53: 1, 62: 1, 64: 2, 65: 1, 67: 1, 88: 1, 90: 1, 92: 1, 95: 1}\n    if len(np.unique(y_true)) > 14:\n        classes.append(99)\n        class_weight[99] = 2\n    y_p = y_preds.reshape(y_true.shape[0], len(classes), order='F')\n\n    # Trasform y_true in dummies\n    y_ohe = pd.get_dummies(y_true)\n    # Normalize rows and limit y_preds to 1e-15, 1-1e-15\n    y_p = np.clip(a=y_p, a_min=1e-15, a_max=1 - 1e-15)\n    # Transform to log\n    y_p_log = np.log(y_p)\n    # Get the log for ones, .values is used to drop the index of DataFrames\n    # Exclude class 99 for now, since there is no class99 in the training set\n    # we gave a special process for that class\n    y_log_ones = np.sum(y_ohe.values * y_p_log, axis=0)\n    # Get the number of positives for each class\n    nb_pos = y_ohe.sum(axis=0).values.astype(float)\n    # Weight average and divide by the number of positives\n    class_arr = np.array([class_weight[k] for k in sorted(class_weight.keys())])\n    y_w = y_log_ones * class_arr / nb_pos\n\n    loss = - np.sum(y_w) / np.sum(class_arr)\n    return 'wloss', loss, False\n\n\ndef save_importances(importances_):\n    mean_gain = importances_[['gain', 'feature']].groupby('feature').mean()\n    importances_['mean_gain'] = importances_['feature'].map(mean_gain['gain'])\n    plt.figure(figsize=(8, 12))\n    sns.barplot(x='gain', y='feature', data=importances_.sort_values('mean_gain', ascending=False))\n    plt.tight_layout()\n    plt.savefig('importances.png')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"e5b519c9151471a624b23b424e74eb27e78bffc0"},"cell_type":"markdown","source":"# 5-fold lgb training"},{"metadata":{"trusted":true,"_uuid":"4a411f6097378ab0c3199064478f4ed4e064b6fb"},"cell_type":"code","source":"\nfolds = StratifiedKFold(n_splits=5, shuffle=True, random_state=1)\nclfs = []\nimportances = pd.DataFrame()\nlgb_params = {\n    'boosting_type': 'gbdt',\n    'objective': 'multiclass',\n    'num_class': 14,\n    'metric': 'multi_logloss',\n    'learning_rate': 0.03,\n    'subsample': .9,\n    'colsample_bytree': 0.5,\n    'reg_alpha': .01,\n    'reg_lambda': .01,\n    'min_split_gain': 0.01,\n    'min_child_weight': 10,\n    'n_estimators': 1000,\n    'silent': -1,\n    'verbose': -1,\n    'max_depth': 3\n}\n\n# Compute weights\nw = y.value_counts()\nweights = {i : np.sum(w) / w[i] for i in w.index}\n\noof_preds = np.zeros((len(full_train), np.unique(y).shape[0]))\nfor fold_, (trn_, val_) in enumerate(folds.split(y, y)):\n    trn_x, trn_y = full_train.iloc[trn_], y.iloc[trn_]\n    val_x, val_y = full_train.iloc[val_], y.iloc[val_]\n\n    clf = lgb.LGBMClassifier(**lgb_params)\n    clf.fit(\n        trn_x, trn_y,\n        eval_set=[(trn_x, trn_y), (val_x, val_y)],\n        eval_metric=lgb_multi_weighted_logloss,\n        verbose=100,\n        early_stopping_rounds=50,\n        sample_weight=trn_y.map(weights)\n    )\n    oof_preds[val_, :] = clf.predict_proba(val_x, num_iteration=clf.best_iteration_)\n    print(multi_weighted_logloss(val_y, oof_preds[val_, :]))\n\n    imp_df = pd.DataFrame()\n    imp_df['feature'] = full_train.columns\n    imp_df['gain'] = clf.feature_importances_\n    imp_df['fold'] = fold_ + 1\n    importances = pd.concat([importances, imp_df], axis=0, sort=False)\n\n    clfs.append(clf)\n\nprint('MULTI WEIGHTED LOG LOSS : %.5f ' % multi_weighted_logloss(y_true=y, y_preds=oof_preds))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5905ad9b53c2d01e2bc33d23feb4b7ccd0e41c50"},"cell_type":"markdown","source":"# Evaluation\nfeatures' importance and confusion matrix"},{"metadata":{"trusted":true,"_uuid":"d8d149b6795271bdce59ee091118bd1681cb3462"},"cell_type":"code","source":"save_importances(importances_=importances)\n# http://scikit-learn.org/stable/modules/generated/sklearn.metrics.confusion_matrix.html\ndef plot_confusion_matrix(cm, classes,\n                          normalize=False,\n                          title='Confusion matrix',\n                          cmap=plt.cm.Blues):\n    \"\"\"\n    This function prints and plots the confusion matrix.\n    Normalization can be applied by setting `normalize=True`.\n    \"\"\"\n    if normalize:\n        cm = cm.astype('float') / cm.sum(axis=1)[:, np.newaxis]\n        print(\"Normalized confusion matrix\")\n    else:\n        print('Confusion matrix, without normalization')\n\n    plt.imshow(cm, interpolation='nearest', cmap=cmap)\n    plt.title(title)\n    plt.colorbar()\n    tick_marks = np.arange(len(classes))\n    plt.xticks(tick_marks, classes, rotation=45)\n    plt.yticks(tick_marks, classes)\n\n    fmt = '.2f' if normalize else 'd'\n    thresh = cm.max() / 2.\n    for i, j in itertools.product(range(cm.shape[0]), range(cm.shape[1])):\n        plt.text(j, i, format(cm[i, j], fmt),\n                 horizontalalignment=\"center\",\n                 color=\"white\" if cm[i, j] > thresh else \"black\")\n\n    plt.ylabel('True label')\n    plt.xlabel('Predicted label')\n    plt.tight_layout()\n    \nunique_y = np.unique(y)\nclass_map = dict()\nfor i,val in enumerate(unique_y):\n    class_map[val] = i\n        \ny_map = np.zeros((y.shape[0],))\ny_map = np.array([class_map[val] for val in y])\n\n# Compute confusion matrix\nfrom sklearn.metrics import confusion_matrix\ncnf_matrix = confusion_matrix(y_map, np.argmax(oof_preds,axis=-1))\nnp.set_printoptions(precision=2)\n\nsample_sub = pd.read_csv('../input/sample_submission.csv')\nclass_names = list(sample_sub.columns[1:-1])\ndel sample_sub;gc.collect()\n\n# Plot non-normalized confusion matrix\nplt.figure(figsize=(12,12))\nfoo = plot_confusion_matrix(cnf_matrix, classes=class_names,normalize=True,\n                      title='Confusion matrix')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"f65492856178f1692de22989481d7dde8afd1e84"},"cell_type":"markdown","source":"# Test Set Predictions"},{"metadata":{"trusted":true,"_uuid":"887b6851259e7e4fddc194f17b13b7110a51c52d"},"cell_type":"code","source":"def predict_chunk(df_, clfs_, meta_, features, train_mean):\n    # Group by object id    \n    agg_ = featurize(df_)\n    # Merge with meta data\n    full_test = agg_.reset_index().merge(\n        right=meta_,\n        how='left',\n        on='object_id'\n    )\n\n    full_test = full_test.fillna(0)\n    # Make predictions\n    preds_ = None\n    for clf in clfs_:\n        if preds_ is None:\n            preds_ = clf.predict_proba(full_test[features]) / len(clfs_)\n        else:\n            preds_ += clf.predict_proba(full_test[features]) / len(clfs_)\n\n    # Compute preds_99 as the proba of class not being any of the others\n    # preds_99 = 0.1 gives 1.769\n    preds_99 = np.ones(preds_.shape[0])\n    for i in range(preds_.shape[1]):\n        preds_99 *= (1 - preds_[:, i])\n\n    # Create DataFrame from predictions\n    preds_df_ = pd.DataFrame(preds_, columns=['class_' + str(s) for s in clfs_[0].classes_])\n    preds_df_['object_id'] = full_test['object_id']\n    preds_df_['class_99'] = 0.14 * preds_99 / np.mean(preds_99) \n    return preds_df_","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"750909d89665158aca660c056e6fe0cd3940f839"},"cell_type":"code","source":"meta_test = pd.read_csv('../input/test_set_metadata.csv')\n# meta_test.set_index('object_id',inplace=True)\nimport time\n\nstart = time.time()\nchunks = 5000000\nremain_df = None\n\nfor i_c, df in enumerate(pd.read_csv('../input/test_set.csv', chunksize=chunks, iterator=True)):\n    # Check object_ids\n    # I believe np.unique keeps the order of group_ids as they appear in the file\n    unique_ids = np.unique(df['object_id'])\n    new_remain_df = df.loc[df['object_id'] == unique_ids[-1]].copy()\n    if remain_df is None:\n        df = df.loc[df['object_id'].isin(unique_ids[:-1])]\n    else:\n        df = pd.concat([remain_df, df.loc[df['object_id'].isin(unique_ids[:-1])]], axis=0)\n    # Create remaining samples df\n    remain_df = new_remain_df\n    preds_df = predict_chunk(df_=df,\n                             clfs_=clfs,\n                             meta_=meta_test,\n                             features=full_train.columns,\n                             train_mean=train_mean)\n\n    if i_c == 0:\n        preds_df.to_csv('predictions.csv', header=True, mode='a', index=False)\n    else:\n        preds_df.to_csv('predictions.csv', header=False, mode='a', index=False)\n\n    del preds_df\n    gc.collect()\n    \n    print('%15d done in %5.1f minutes' % (chunks * (i_c + 1), (time.time() - start) / 60), flush=True)\n\n# Compute last object in remain_df\npreds_df = predict_chunk(df_=remain_df,\n                         clfs_=clfs,\n                         meta_=meta_test,\n                         features=full_train.columns,\n                         train_mean=train_mean)\n\npreds_df.to_csv('predictions.csv', header=False, mode='a', index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"71e0c64350811ff7b0c34af9bb85bdc3f3c51068"},"cell_type":"code","source":"z = pd.read_csv('predictions.csv')\nprint(\"Shape BEFORE grouping:\",z.shape)\nz = z.groupby('object_id').mean()\nprint(\"Shape AFTER grouping:\",z.shape)\nz.to_csv('single_predictions.csv', index=True)","execution_count":null,"outputs":[]}],"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}