""" This script is forked from iprapas's notebook https://www.kaggle.com/iprapas/ideas-from-kernels-and-discussion-lb-1-135 # https://www.kaggle.com/ogrellier/plasticc-in-a-kernel-meta-and-data # https://www.kaggle.com/c/PLAsTiCC-2018/discussion/70908 # https://www.kaggle.com/meaninglesslives/simple-neural-net-for-time-series-classification # """ import sys, os import argparse import time from datetime import datetime as dt import gc; gc.enable() from functools import partial, wraps import numpy as np # linear algebra import pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv) np.warnings.filterwarnings('ignore') from sklearn.model_selection import StratifiedKFold from tsfresh.feature_extraction import extract_features from xgboost import XGBClassifier from lightgbm import LGBMClassifier def featurize(df, fcp, aggs, n_jobs=4): """ Extracting Features from train set Features from olivier's kernel very smart and powerful feature that is generously given here https://www.kaggle.com/c/PLAsTiCC-2018/discussion/69696#410538 per passband features with tsfresh library. fft features added to capture periodicity https://www.kaggle.com/c/PLAsTiCC-2018/discussion/70346#415506 """ df['flux_ratio_sq'] = np.power(df['flux'].values / df['flux_err'].values, 2.0) df['flux_by_flux_ratio_sq'] = df['flux'].values * df['flux_ratio_sq'].values # train[detected==1, mjd_diff:=max(mjd)-min(mjd), by=object_id] agg_df = df.groupby('object_id').agg(aggs) agg_df.columns = [ '{}_{}'.format(k, agg) for k in aggs.keys() for agg in aggs[k]] agg_df['flux_w_mean'] = agg_df['flux_by_flux_ratio_sq_sum'].values / agg_df['flux_ratio_sq_sum'].values agg_df['flux_diff'] = agg_df['flux_max'].values - agg_df['flux_min'].values agg_df['flux_dif2'] = agg_df['flux_diff'].values / agg_df['flux_mean'].values agg_df['flux_dif3'] = agg_df['flux_diff'].values / agg_df['flux_w_mean'].values # Add more features with 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=n_jobs) # Add smart feature that is suggested here https://www.kaggle.com/c/PLAsTiCC-2018/discussion/69696#410538 # dt[detected==1, mjd_diff:=max(mjd)-min(mjd), by=object_id] df_det = df[df['detected']==1].copy() agg_df_mjd = extract_features(df_det, column_id='object_id', column_value='mjd', default_fc_parameters={'maximum':None, 'minimum':None}, n_jobs=n_jobs) agg_df_mjd['mjd_diff_det'] = agg_df_mjd['mjd__maximum'].values - agg_df_mjd['mjd__minimum'].values del agg_df_mjd['mjd__maximum'], agg_df_mjd['mjd__minimum'] agg_df_ts = pd.merge(agg_df_ts, agg_df_mjd, on='id') # tsfresh returns a dataframe with an index name='id' agg_df_ts.index.rename('object_id', inplace=True) agg_df = pd.merge(agg_df, agg_df_ts, on='object_id') return agg_df def multi_weighted_logloss(y_true, y_preds, classes, class_weights): """ refactor from @author olivier https://www.kaggle.com/ogrellier multi logloss for PLAsTiCC challenge """ y_p = y_preds.reshape(y_true.shape[0], len(classes), order='F') # Trasform y_true in dummies y_ohe = pd.get_dummies(y_true) # Normalize rows and limit y_preds to 1e-15, 1-1e-15 y_p = np.clip(a=y_p, a_min=1e-15, a_max=1 - 1e-15) # Transform to log y_p_log = np.log(y_p) # Get the log for ones, .values is used to drop the index of DataFrames # Exclude class 99 for now, since there is no class99 in the training set # we gave a special process for that class y_log_ones = np.sum(y_ohe.values * y_p_log, axis=0) # Get the number of positives for each class nb_pos = y_ohe.sum(axis=0).values.astype(float) # Weight average and divide by the number of positives class_arr = np.array([class_weights[k] for k in sorted(class_weights.keys())]) y_w = y_log_ones * class_arr / nb_pos loss = - np.sum(y_w) / np.sum(class_arr) return loss def lgbm_multi_weighted_logloss(y_true, y_preds): """ refactor from @author olivier https://www.kaggle.com/ogrellier multi logloss for PLAsTiCC challenge """ # Taken from Giba's topic : https://www.kaggle.com/titericz # https://www.kaggle.com/c/PLAsTiCC-2018/discussion/67194 # with Kyle Boone's post https://www.kaggle.com/kyleboone classes = [6, 15, 16, 42, 52, 53, 62, 64, 65, 67, 88, 90, 92, 95] class_weights = {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} loss = multi_weighted_logloss(y_true, y_preds, classes, class_weights) return 'wloss', loss, False def save_importances(importances_): mean_gain = importances_[['gain', 'feature']].groupby('feature').mean() importances_['mean_gain'] = importances_['feature'].map(mean_gain['gain']) return importances_ def lgbm_modeling_cross_validation(params, full_train, y, classes, class_weights, nr_fold=5, random_state=1): # Compute weights w = y.value_counts() weights = {i : np.sum(w) / w[i] for i in w.index} clfs = [] importances = pd.DataFrame() folds = StratifiedKFold(n_splits=nr_fold, shuffle=True, random_state=random_state) oof_preds = np.zeros((len(full_train), np.unique(y).shape[0])) for fold_, (trn_, val_) in enumerate(folds.split(y, y)): trn_x, trn_y = full_train.iloc[trn_], y.iloc[trn_] val_x, val_y = full_train.iloc[val_], y.iloc[val_] clf = LGBMClassifier(**params) clf.fit( trn_x, trn_y, eval_set=[(trn_x, trn_y), (val_x, val_y)], eval_metric=lgbm_multi_weighted_logloss, verbose=100, early_stopping_rounds=50, sample_weight=trn_y.map(weights) ) clfs.append(clf) oof_preds[val_, :] = clf.predict_proba(val_x, num_iteration=clf.best_iteration_) print('no {}-fold loss: {}'.format(fold_ + 1, multi_weighted_logloss(val_y, oof_preds[val_, :], classes, class_weights))) imp_df = pd.DataFrame({ 'feature': full_train.columns, 'gain': clf.feature_importances_, 'fold': [fold_ + 1] * len(full_train.columns), }) importances = pd.concat([importances, imp_df], axis=0, sort=False) score = multi_weighted_logloss(y_true=y, y_preds=oof_preds, classes=classes, class_weights=class_weights) print('MULTI WEIGHTED LOG LOSS: {:.5f}'.format(score)) df_importances = save_importances(importances_=importances) df_importances.to_csv('xgb_importances.csv', index=False) return clfs, score def predict_chunk(df_, clfs_, meta_, features, featurize_configs, train_mean): # Group by object id agg_ = featurize(df_, featurize_configs['fcp'], featurize_configs['aggs']) # Merge with meta data agg_.reset_index(inplace=True) full_test = agg_.merge(right=meta_, how='left', on='object_id') full_test.fillna(0, inplace=True) # Make predictions preds_ = None for clf in clfs_: if preds_ is None: preds_ = clf.predict_proba(full_test[features]) / len(clfs_) else: preds_ += clf.predict_proba(full_test[features]) / len(clfs_) # Compute preds_99 as the proba of class not being any of the others # preds_99 = 0.1 gives 1.769 preds_99 = np.ones(preds_.shape[0]) for i in range(preds_.shape[1]): preds_99 *= (1 - preds_[:, i]) # Create DataFrame from predictions preds_df_ = pd.DataFrame(preds_, columns=['class_{}'.format(s) for s in clfs_[0].classes_]) preds_df_['object_id'] = full_test['object_id'] preds_df_['class_99'] = 0.14 * preds_99 / np.mean(preds_99) return preds_df_ def process_test(clfs, features, featurize_configs, train_mean, filename='predictions.csv', chunks=5000000): start = time.time() meta_test = pd.read_csv('../input/test_set_metadata.csv') # meta_test.set_index('object_id',inplace=True) remain_df = None for i_c, df in enumerate(pd.read_csv('../input/test_set.csv', chunksize=chunks, iterator=True)): # Check object_ids # I believe np.unique keeps the order of group_ids as they appear in the file unique_ids = np.unique(df['object_id']) new_remain_df = df.loc[df['object_id'] == unique_ids[-1]].copy() if remain_df is None: df = df.loc[df['object_id'].isin(unique_ids[:-1])] else: df = pd.concat([remain_df, df.loc[df['object_id'].isin(unique_ids[:-1])]], axis=0) # Create remaining samples df remain_df = new_remain_df preds_df = predict_chunk(df_=df, clfs_=clfs, meta_=meta_test, features=features, featurize_configs=featurize_configs, train_mean=train_mean) if i_c == 0: preds_df.to_csv(filename, header=True, mode='a', index=False) else: preds_df.to_csv(filename, header=False, mode='a', index=False) del preds_df gc.collect() print('{:15d} done in {:5.1f} minutes' .format( chunks * (i_c + 1), (time.time() - start) / 60), flush=True) # Compute last object in remain_df preds_df = predict_chunk(df_=remain_df, clfs_=clfs, meta_=meta_test, features=features, featurize_configs=featurize_configs, train_mean=train_mean) preds_df.to_csv(filename, header=False, mode='a', index=False) return def main(argc, argv): # Features to compute with tsfresh library. Fft coefficient is meant to capture periodicity fcp = {'fft_coefficient': [ {'coeff': 0, 'attr': 'abs'}, {'coeff': 1, 'attr': 'abs'}], 'kurtosis' : None, 'skewness' : None} aggs = { 'flux': ['min', 'max', 'mean', 'median', 'std', 'skew'], 'flux_err': ['min', 'max', 'mean', 'median', 'std', 'skew'], 'detected': ['mean'], 'flux_ratio_sq':['sum','skew'], 'flux_by_flux_ratio_sq':['sum','skew'], } best_params = { 'device': 'cpu', 'objective': 'multiclass', 'num_class': 14, 'boosting_type': 'gbdt', 'n_jobs': -1, 'max_depth': 6, 'n_estimators': 1000, 'subsample_freq': 2, 'subsample_for_bin': 5000, 'min_data_per_group': 100, 'max_cat_to_onehot': 4, 'cat_l2': 1.125, 'cat_smooth': 100.0, 'max_cat_threshold': 32, 'metric_freq': 10, 'verbosity': -1, 'metric': 'multi_logloss', 'xgboost_dart_mode': False, 'uniform_drop': False, 'colsample_bytree': 0.65, 'drop_rate': 0.25, 'learning_rate': 0.03, 'max_drop': 25, 'min_child_samples': 10, 'min_child_weight': 100.0, 'min_split_gain': 0.0008, 'num_leaves': 11, 'reg_alpha': 0.1, 'reg_lambda': 1.0, 'skip_drop': 0.75, 'subsample': 0.6} train = pd.read_csv('../input/training_set.csv') agg_train = featurize(train, fcp, aggs) meta_train = pd.read_csv('../input/training_set_metadata.csv') meta_train.head() full_train = agg_train.reset_index().merge(right=meta_train, how='outer', on='object_id' ) if 'target' in full_train: y = full_train['target'] del full_train['target'] classes = sorted(y.unique()) # Taken from Giba's topic : https://www.kaggle.com/titericz # https://www.kaggle.com/c/PLAsTiCC-2018/discussion/67194 # with Kyle Boone's post https://www.kaggle.com/kyleboone class_weights = {c: 1 for c in classes} class_weights.update({c:2 for c in [64, 15]}) print('Unique classes : {}, {}'.format(len(classes), classes)) print(class_weights) #sanity check: classes = [6, 15, 16, 42, 52, 53, 62, 64, 65, 67, 88, 90, 92, 95] #sanity check: class_weights = {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} #if len(np.unique(y_true)) > 14: # classes.append(99) # class_weights[99] = 2 if 'object_id' in full_train: oof_df = full_train[['object_id']] del full_train['object_id'], full_train['distmod'], full_train['hostgal_specz'] del full_train['ra'], full_train['decl'], full_train['gal_l'], full_train['gal_b'], full_train['ddf'] train_mean = full_train.mean(axis=0) full_train.fillna(0, inplace=True) # modeling from CV eval_func = partial(lgbm_modeling_cross_validation, full_train=full_train, y=y, classes=classes, class_weights=class_weights, nr_fold=5, random_state=1) clfs, score = eval_func(best_params) filename = 'subm_{:.6f}_{}.csv'.format(score, dt.now().strftime('%Y-%m-%d-%H-%M')) print('save to {}'.format(filename)) # TEST process_test(clfs, features=full_train.columns, featurize_configs={'fcp': fcp, 'aggs': aggs}, train_mean=train_mean, filename=filename, chunks=5000000) z = pd.read_csv(filename) print("Shape BEFORE grouping: {}".format(z.shape)) z = z.groupby('object_id').mean() print("Shape AFTER grouping: {}".format(z.shape)) z.to_csv('single_{}'.format(filename), index=True) if __name__ == '__main__': main(len(sys.argv), sys.argv)