{"cells":[{"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\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 read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        #print(os.path.join(dirname, filename))\n        continue\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session\n\n# Imports\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nplt.style.use('seaborn-white')\n%matplotlib inline\n\nfrom sklearn.model_selection import train_test_split, GridSearchCV\nimport lightgbm as lgbm\nfrom lightgbm import LGBMRegressor\nfrom xgboost import XGBRegressor","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* **Version 1:** Usage of entire data with no special feature selection, model is LGBM with ~default parameters.\n* **Version 2:** Version 1 plus Grid Search of best LGBM parameters."},{"metadata":{},"cell_type":"markdown","source":"# Data"},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"PATH = '/kaggle/input/predict-volcanic-eruptions-ingv-oe/'\n\n!echo 'Files: train and test'\n!ls -l /kaggle/input/predict-volcanic-eruptions-ingv-oe/train/ | wc -l\n!ls -l /kaggle/input/predict-volcanic-eruptions-ingv-oe/test/  | wc -l\n\ntrain_files = []\ntest_files  = []\n\nfor file in os.listdir(PATH+'/train/'):\n    train_files.append(file)\n    \nfor file in os.listdir(PATH+'/test/'):\n    test_files.append(file)\n    \nprint('Number of train files: {}'.format(len(train_files)))\nprint('Number of test  files: {}'.format(len(test_files )))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# Train file\ntrain = pd.read_csv(PATH+'train.csv')\n\n# Submission file\nsample_submission = pd.read_csv(PATH+'sample_submission.csv')\n\ntest_files = []\nfor dirname, _, filenames in os.walk(PATH+'/test/'):\n    for filename in filenames:\n        test_files.append(filename[:-4]) # without .csv extension\n        \ntest = pd.DataFrame(test_files, columns=[\"segment_id\"])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## EDA\n### Time to eruption distribution  \n* Looks almost uniform"},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(train['time_to_eruption'], \n             hist=True, \n             kde=True, \n             bins=100, \n             color = 'blue', \n             hist_kws={'edgecolor':'black'})","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# Check min and max in segment_id by sorting\n# train.sort_values('time_to_eruption', axis=0, ascending=True)\ndisplay(train.sort_values('time_to_eruption', axis=0, ascending=True).iloc[[0,-1],:])\n\nsegment_id_min =  601524801\nsegment_id_max = 1923243961\n\ndf_segment_id_min = pd.read_csv(PATH+'/train/'+str(segment_id_min)+'.csv')\ndf_segment_id_max = pd.read_csv(PATH+'/train/'+str(segment_id_max)+'.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df_segment_id_min.head(3)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Sensors"},{"metadata":{"trusted":true},"cell_type":"code","source":"df_segment_id_min.plot(figsize=(20,20),\n                       subplots=True, \n                       layout=(10,1),\n                       rot=0, \n                       lw=1, \n                       #colormap='jet',\n                       title='segment_id #601524801 (min)'\n                      )\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df_segment_id_max.plot(figsize=(20,20),\n                       subplots=True, \n                       layout=(10,1),\n                       rot=0, \n                       lw=2, \n                       #colormap='jet',\n                       title='segment_id #1923243961 (max)'\n                      )\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"So,\n  * sensors 2, 3, 8 are out (empty) in `id_min`; some spikes suggest something abnormal\n  * sensors 1 and 5   are out (empty) in `id_max`; no spikes, smoothness suggests something normal"},{"metadata":{},"cell_type":"markdown","source":"## Features\n* So, as we see, we able to get outputs for sensor 1-10 from each {segment_id}.csv.  \n* Let's build features on that. With a lot of notebooks available we can select some findings there (as below).  \n* Please upvote [this](https://www.kaggle.com/isaienkov/ingv-volcanic-eruption-prediction-eda-modeling) amazing notebook from [Kostiantyn Isaienkov](https://www.kaggle.com/isaienkov) where this clear function below is taken from. "},{"metadata":{"trusted":true},"cell_type":"code","source":"def build_features(signal, ts, sensor_id):\n    X = pd.DataFrame()\n    f = np.fft.fft(signal)\n    f_real = np.real(f)\n    X.loc[ts, f'{sensor_id}_sum']       = signal.sum()\n    X.loc[ts, f'{sensor_id}_mean']      = signal.mean()\n    X.loc[ts, f'{sensor_id}_std']       = signal.std()\n    X.loc[ts, f'{sensor_id}_var']       = signal.var() \n    X.loc[ts, f'{sensor_id}_max']       = signal.max()\n    X.loc[ts, f'{sensor_id}_min']       = signal.min()\n    X.loc[ts, f'{sensor_id}_skew']      = signal.skew()\n    X.loc[ts, f'{sensor_id}_mad']       = signal.mad()\n    X.loc[ts, f'{sensor_id}_kurtosis']  = signal.kurtosis()\n    X.loc[ts, f'{sensor_id}_quantile99']= np.quantile(signal, 0.99)\n    X.loc[ts, f'{sensor_id}_quantile95']= np.quantile(signal, 0.95)\n    X.loc[ts, f'{sensor_id}_quantile85']= np.quantile(signal, 0.85)\n    X.loc[ts, f'{sensor_id}_quantile75']= np.quantile(signal, 0.75)\n    X.loc[ts, f'{sensor_id}_quantile55']= np.quantile(signal, 0.55)\n    X.loc[ts, f'{sensor_id}_quantile45']= np.quantile(signal, 0.45) \n    X.loc[ts, f'{sensor_id}_quantile25']= np.quantile(signal, 0.25) \n    X.loc[ts, f'{sensor_id}_quantile15']= np.quantile(signal, 0.15) \n    X.loc[ts, f'{sensor_id}_quantile05']= np.quantile(signal, 0.05)\n    X.loc[ts, f'{sensor_id}_quantile01']= np.quantile(signal, 0.01)\n    X.loc[ts, f'{sensor_id}_fft_real_mean']= f_real.mean()\n    X.loc[ts, f'{sensor_id}_fft_real_std'] = f_real.std()\n    X.loc[ts, f'{sensor_id}_fft_real_max'] = f_real.max()\n    X.loc[ts, f'{sensor_id}_fft_real_min'] = f_real.min()\n\n    return X","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We will perform some statistics on each column (sensors 1-10) obtained from {segment_id}.csv files."},{"metadata":{},"cell_type":"markdown","source":"## Dataset"},{"metadata":{},"cell_type":"markdown","source":"### Train"},{"metadata":{"trusted":true},"cell_type":"code","source":"train_set = list()\nseg=0\n\nfor seg, segment_id in enumerate(train.segment_id):\n    signals = pd.read_csv(PATH+'/train/'+str(segment_id)+'.csv')\n    train_row = []\n    \n    if seg%200 == 0:\n        print('Processing segment_id={}'.format(seg))\n        \n    for sensor in range(0, 10):\n        sensor_id = f'sensor_{sensor+1}'\n        train_row.append(build_features(signals[sensor_id].fillna(0), segment_id, sensor_id))\n        \n    train_row = pd.concat(train_row, axis=1)\n    train_set.append(train_row)\n    seg+=1\n    \ntrain_set = pd.concat(train_set)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_set = train_set.reset_index()                                        \ntrain_set = train_set.rename(columns={'index': 'segment_id'}) # change column index->segment_id\ntrain_set = pd.merge(train_set, train, on='segment_id')       # merge with our train by segment_id column","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_set.head(3)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Test"},{"metadata":{"trusted":true},"cell_type":"code","source":"test_set = list()\nseg=0\n\nfor seg, segment_id in enumerate(test.segment_id):\n    signals = pd.read_csv(PATH+'/test/'+str(segment_id)+'.csv')\n    test_row = []\n    \n    if seg%200 == 0:\n        print('Processing segment_id={}'.format(seg))\n        \n    for sensor in range(0, 10):\n        sensor_id = f'sensor_{sensor+1}'\n        test_row.append(build_features(signals[sensor_id].fillna(0), segment_id, sensor_id))\n        \n    test_row = pd.concat(test_row, axis=1)\n    test_set.append(test_row)\n    seg+=1\n    \ntest_set = pd.concat(test_set)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_set = test_set.reset_index()                                        \ntest_set = test_set.rename(columns={'index': 'segment_id'}) # change column index->segment_id\ntest_set = pd.merge(test_set, test, on='segment_id')        # merge with our train by segment_id column","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_set.head(3)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Train/Test split"},{"metadata":{"trusted":true},"cell_type":"code","source":"X = train_set.drop(['segment_id', 'time_to_eruption'], axis=1)\ny = train_set['time_to_eruption']\n\nX_train, X_valid, y_train, y_valid = train_test_split(X, y, \n                                                      test_size=0.2, \n                                                      random_state=42)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Models\n\nWe will use:\n* LGBM\n* XGBoost"},{"metadata":{},"cell_type":"markdown","source":"## LGBM\n\nLight GBM is a gradient boosting framework that uses tree based learning algorithm.\n\nSome parameters:  \n* **boosting_type** – ‘gbdt’ (default), ‘dart’, ‘goss’, ‘rf’. Default gbdt is gradient boosting decision tree.\n* **num_leaves** – Maximum tree leaves (default=31) for base learners.\n* **max_depth** – Maximum tree depth (default=-1) for base learners, <=0 means no limit.\n* **learning_rate**  – Boosting learning rate (default=0.1). \n* **n_estimators** – Number of boosted trees (default=100) to fit.\n* **subsample_for_bin** – Number of samples for constructing bins (default=200000).\n* **objective** – Specify the learning task and the corresponding learning objective, e.g, ‘binary’ or ‘multiclass’.\n* **min_split_gain** – Minimum loss reduction (default=0.) required to make a further partition on a leaf node of tree.\n* **min_child_weight** – Minimum sum of instance weight (hessian) needed in a child/leaf (default=1e-3).\n* **min_child_samples** – Minimum number of data needed in a child/leaf (default=20).\n* **subsample** – Subsample ratio of the training instance (default=1.).\n* **subsample_freq** – Frequence of subsample (default=0), <=0 means no enable.\n* **colsample_bytree** – Subsample ratio (default=1.) of columns when constructing each tree.\n* **reg_alpha** – L1 regularization term on weights (default=0.).\n* **reg_lambda** – L2 regularization term on weights (default=0.).\n* **random_state** – Random number seed. \n* **n_jobs** – Number of parallel threads (default=-1).\n* **silent** – Whether to print messages while running boosting (default=True).  \n\n\n* **min_data_in_leaf** – Minimal number of data in one leaf (default=20). Can be used to deal with over-fitting\n* **feature_fraction** – LightGBM will randomly select part of features (default=1.0) on each iteration (tree) if feature_fraction smaller than 1.0. For example, if you set it to 0.8, LightGBM will select 80% of features before training each tree. Can be used to speed up training and to deal with over-fitting.\n* **bagging_fraction** – LightGBM will randomly select part of data (default=1.0) without resampling. Can be used to speed up training and to deal with over-fitting\n\n        https://lightgbm.readthedocs.io/en/latest/pythonapi/lightgbm.LGBMRegressor.html\n        https://lightgbm.readthedocs.io/en/latest/Parameters.html"},{"metadata":{},"cell_type":"markdown","source":"### Parameters"},{"metadata":{"trusted":true},"cell_type":"code","source":"# Default parameters\nparams = {\n    'boosting_type': 'gbdt', \n    'num_leaves': 31,\n    'max_depth': -1,\n    'learning_rate': 0.1,\n    'n_estimators': 100, \n    'subsample_for_bin': 200, # 200000 is default \n    #'objective': 'binary'\n    'min_split_gain': 0.5,    # 0.0 is default \n    'min_child_weight': 1e-3, \n    'min_child_samples': 20,\n    'subsample': 1,\n    'colsample_bytree': 1.0,\n    'min_data_in_leaf': 20,\n    'feature_fraction': 1.0,\n    'bagging_fraction': 1.0,\n    'random_state': 42\n    #'device': 'cpu', # you can use GPU to achieve faster learning\n}\n        \n# Initiate classifier to use\nmodel_lgbm_regr = LGBMRegressor(boosting_type = params['boosting_type'], \n                                num_leaves = params['num_leaves'],\n                                max_depth = params['max_depth'],\n                                learning_rate = params['learning_rate'],\n                                n_estimators = params['n_estimators'],\n                                subsample_for_bin = params['subsample_for_bin'],\n                                #objective = params['objective'],\n                                min_split_gain = params['min_split_gain'], \n                                min_child_weight = params['min_child_weight'], \n                                min_child_samples = params['min_child_samples'],\n                                subsample = params['subsample'],\n                                colsample_bytree = params['colsample_bytree'],\n                                min_data_in_leaf = params['min_data_in_leaf'],\n                                feature_fraction = params['feature_fraction'],\n                                bagging_fraction = params['bagging_fraction'],\n                                random_state = params['random_state'],\n                                #n_jobs = 5, \n                                silent = True\n                               )\n\n# To view the default model parameters:\nmodel_lgbm_regr.get_params().keys()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Grid Search"},{"metadata":{"trusted":true},"cell_type":"code","source":"gridParams = {\n    'boosting_type' : ['gbdt', 'dart'], # for better accuracy might try dart, check both\n    'num_leaves': [10,31],        # large num_leaves helps improve accuracy but might lead to over-fitting\n    'max_depth': [30,10,-1],\n    'learning_rate': [0.01, 0.05, 0.1],\n    'n_estimators': [50, 100, 200],\n    'subsample_for_bin': [200,500],\n    'objective' : ['binary'],\n    'min_split_gain': [0.5,0.8],    # 0.0 is default \n    'min_child_weight': [1e-3,1e-1,1e-2], \n    'min_child_samples': [20],\n    'subsample': [0.6, 0.8, 1.0],\n    'colsample_bytree': [0.6, 0.8, 1.0],\n    'min_data_in_leaf': [5, 10, 20],\n    'feature_fraction': [0.8, 1.0],\n    'bagging_fraction': [0.8, 1.0],\n    'random_state' : [42],\n    }\n\ngrid = GridSearchCV(model_lgbm_regr, \n                    gridParams, \n                    cv=5, \n                    verbose=1, \n                    n_jobs=-1)\n\n# Run the grid\ngrid.fit(X_train, y_train)\n\n# Print the best parameters found\nprint(grid.best_params_)\nprint(grid.best_score_)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Best model selection"},{"metadata":{"trusted":true},"cell_type":"code","source":"params['boosting_type'] = grid.best_params_['boosting_type']\nparams['num_leaves'] = grid.best_params_['num_leaves']\nparams['max_depth'] = grid.best_params_['max_depth']\nparams['learning_rate'] = grid.best_params_['learning_rate']\nparams['n_estimators'] = grid.best_params_['n_estimators']\nparams['subsample_for_bin'] = grid.best_params_['subsample_for_bin']\nparams['objective'] = grid.best_params_['objective']\nparams['min_split_gain'] = grid.best_params_['min_split_gain']\nparams['min_child_weight'] = grid.best_params_['min_child_weight']\nparams['min_child_samples'] = grid.best_params_['min_child_samples']\nparams['subsample'] = grid.best_params_['subsample']\nparams['colsample_bytree'] = grid.best_params_['colsample_bytree']\nparams['min_data_in_leaf'] = grid.best_params_['min_data_in_leaf']\nparams['feature_fraction'] = grid.best_params_['feature_fraction']\nparams['bagging_fraction'] = grid.best_params_['bagging_fraction']\nparams['random_state'] = grid.best_params_['random_state']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#model_lgbm_regr.fit(X_train, y_train, \n#                    eval_set= [(X_train, y_train), (X_valid, y_valid)], \n#                    eval_metric=\"mae\", \n#                    verbose=200, \n#                    early_stopping_rounds=50\n#                   )","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dtrain = lgbm.Dataset(X_train, label=y_train)\ndvalid = lgbm.Dataset(X_valid, label=y_valid) ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model_lgbm = lgbm.train(params, \n                        train_set=dtrain, \n                        num_boost_round=100, \n                        valid_sets=[dvalid, dtrain], \n                        early_stopping_rounds=20, \n                        verbose_eval=4\n                       )","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model_lgbm","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Plot importance"},{"metadata":{"trusted":true},"cell_type":"code","source":"#lgbm.plot_importance(model_lgbm)\n#plt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Predict"},{"metadata":{"trusted":true},"cell_type":"code","source":"#predictions = model_lgbm_regr.predict(test_set.drop(columns=['segment_id']))\npredictions  = model_lgbm.predict(test_set.drop(columns=['segment_id']))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Submission"},{"metadata":{"trusted":true},"cell_type":"code","source":"sample_submission","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission = pd.DataFrame()\nsubmission['segment_id'] = test_set[\"segment_id\"]\nsubmission['time_to_eruption'] = predictions\nsubmission.to_csv('submission.csv', header=True, index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"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":4,"nbformat_minor":4}