{"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_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"<h1 align=\"center\"> LANL Earthquake Prediction</h1>","metadata":{}},{"cell_type":"markdown","source":"<img style=\"margin: 0 auto;\" src=\"https://storage.googleapis.com/kaggle-media/competitions/LANL/nik-shuliahin-585307-unsplash.jpg\" width=\"450\" height=\"200\" style=\"align:center\"/>","metadata":{}},{"cell_type":"markdown","source":"# Table of Contents\n[1. Introduction](#sec1) <br>\n[2. Import Libraries](#sec2)<br>\n[3. Data Processing](#sec3)<br>\n[4. Evaluating Some Models](#sec4)<br>\n[5. Hyperparameter Tuning](#sec5)<br>\n[6. Building Final Models](#sec6)<br>\n[7. Prediction](#sec7)\n","metadata":{}},{"cell_type":"markdown","source":"# 1. Introduction\n<a id='sec1'></a>\nThis competition aims to use seismic signals to predict the timing of laboratory earthquakes. The data is from an experimental set-up. The acoustic_data input signal is used to predict the time remaining before the next laboratory earthquake (time_to_failure).<br>\nThe training data is a single continuous segment of experimental data. The test folder contains many small segments, which do not represent a continuous segment of the experiment. Predictions cannot be assumed to follow the same regular pattern seen in the training file.","metadata":{}},{"cell_type":"markdown","source":"# 2. Import Libraries\n<a id='sec2'></a>","metadata":{}},{"cell_type":"code","source":"import os\nimport glob\nimport numpy as np\nfrom random import seed\nfrom random import randint\n\nfrom tqdm import tqdm\n# Fix seeds\nfrom numpy.random import seed\nseed(639)\n\nfrom scipy.signal import hilbert\nfrom scipy.signal import convolve\nfrom scipy.signal import hann\n\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.model_selection import KFold, cross_val_score\n\nfrom sklearn.linear_model import LinearRegression, Ridge, Lasso\nfrom sklearn.neighbors import KNeighborsRegressor\nfrom sklearn.tree import DecisionTreeRegressor\nfrom sklearn.ensemble import (RandomForestRegressor, GradientBoostingRegressor, \n                              AdaBoostRegressor)\n#from sklearn.ensemble import AdaBoostRegressor\n\nfrom sklearn.ensemble import VotingRegressor\n\nfrom sklearn.model_selection import RandomizedSearchCV\n\nfrom sklearn.model_selection import RepeatedKFold\nfrom sklearn.model_selection import GridSearchCV\nimport seaborn as sns\n\n\nimport pandas as pd\n\nimport matplotlib.pyplot as plt\nfrom tqdm import tqdm_notebook\nimport warnings\nwarnings.filterwarnings('ignore')\npd.options.display.precision = 15","metadata":{"execution":{"iopub.status.busy":"2021-07-03T19:48:32.572610Z","iopub.execute_input":"2021-07-03T19:48:32.573196Z","iopub.status.idle":"2021-07-03T19:48:34.172463Z","shell.execute_reply.started":"2021-07-03T19:48:32.573108Z","shell.execute_reply":"2021-07-03T19:48:34.171235Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 3. Data Processing\n<a id='sec3'></a>","metadata":{}},{"cell_type":"code","source":"# Load training data\nINPATH = '../input/LANL-Earthquake-Prediction/'\ntrain_df = pd.read_csv(os.path.join( INPATH, 'train.csv'),\\\n                       dtype={'acoustic_data': np.int16, 'time_to_failure': np.float32})\n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T19:48:34.174568Z","iopub.execute_input":"2021-07-03T19:48:34.174951Z","iopub.status.idle":"2021-07-03T19:52:07.697146Z","shell.execute_reply.started":"2021-07-03T19:48:34.174916Z","shell.execute_reply":"2021-07-03T19:52:07.694821Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### 3.1. Useful Functions\nWe'll use the following functions for data processing and feature engineering.","metadata":{}},{"cell_type":"code","source":"# sta/lta function from obspy: https://docs.obspy.org/_modules/obspy/signal/trigger.html\ndef recursive_sta_lta(a, nsta, nlta):\n    '''\n    a: seismic trace\n    nsta: Length of short time average window in samples\n    nlta: Length of long time average window in samples\n    '''\n    try:\n        a = a.tolist()\n    except Exception:\n        pass\n    ndat = len(a)\n    # compute the short time average (STA) and long time average (LTA)\n    # given by Evans and Allen\n    csta = 1. / nsta\n    clta = 1. / nlta\n    sta = 0.\n    lta = 1e-99  # avoid zero division\n    charfct = [0.0] * len(a)\n    icsta = 1 - csta\n    iclta = 1 - clta\n    for i in range(1, ndat):\n        sq = a[i] ** 2\n        sta = csta * sq + icsta * sta\n        lta = clta * sq + iclta * lta\n        charfct[i] = sta / lta\n        if i < nlta:\n            charfct[i] = 0.\n    return np.array(charfct)","metadata":{"execution":{"iopub.status.busy":"2021-07-03T19:52:07.702560Z","iopub.execute_input":"2021-07-03T19:52:07.703344Z","iopub.status.idle":"2021-07-03T19:52:07.716401Z","shell.execute_reply.started":"2021-07-03T19:52:07.703288Z","shell.execute_reply":"2021-07-03T19:52:07.715081Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def add_linear_trend(arr):\n    idx = np.array(range(len(arr)))\n    lr = LinearRegression()\n    lr.fit(idx.reshape(-1, 1), arr)\n    return lr.coef_[0]","metadata":{"execution":{"iopub.status.busy":"2021-07-03T19:52:07.717791Z","iopub.execute_input":"2021-07-03T19:52:07.718083Z","iopub.status.idle":"2021-07-03T19:52:07.737084Z","shell.execute_reply.started":"2021-07-03T19:52:07.718055Z","shell.execute_reply":"2021-07-03T19:52:07.735765Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#https://www.kaggle.com/gpreda/lanl-earthquake-eda-and-prediction\ndef create_features(seg_id, seg, X):\n    xc = pd.Series(seg['acoustic_data'].values)\n    zc = np.fft.fft(xc)\n    \n    X.loc[seg_id, 'mean'] = xc.mean()\n    X.loc[seg_id, 'std'] = xc.std()\n    X.loc[seg_id, 'max'] = xc.max()\n    X.loc[seg_id, 'min'] = xc.min()\n    \n    X.loc[seg_id, 'mad'] = xc.mad()\n    X.loc[seg_id, 'kurt'] = xc.kurtosis()\n    X.loc[seg_id, 'skew'] = xc.skew()\n    X.loc[seg_id, 'med'] = xc.median()\n    \n    #FFT transform values\n    realFFT = np.real(zc)\n    imagFFT = np.imag(zc)\n    X.loc[seg_id, 'Rmean'] = realFFT.mean()\n    X.loc[seg_id, 'Rstd'] = realFFT.std()\n    X.loc[seg_id, 'Rmax'] = realFFT.max()\n    X.loc[seg_id, 'Rmin'] = realFFT.min()\n    X.loc[seg_id, 'Imean'] = imagFFT.mean()\n    X.loc[seg_id, 'Istd'] = imagFFT.std()\n    X.loc[seg_id, 'Imax'] = imagFFT.max()\n    X.loc[seg_id, 'Imin'] = imagFFT.min()\n    \n    X.loc[seg_id, 'mean_change_abs'] = np.mean(np.diff(xc))\n       \n    X.loc[seg_id, 'abs_max'] = np.abs(xc).max()\n    X.loc[seg_id, 'abs_min'] = np.abs(xc).min()\n    X.loc[seg_id, 'abs_std'] = np.abs(xc).std()\n    X.loc[seg_id, 'abs_mean'] = np.abs(xc).mean()\n    \n    for value in [10000, 50000]:\n        #value = randint(10000, 150000)\n        X.loc[seg_id, 'std_first_%s' %value] = xc[:value].std()\n        X.loc[seg_id, 'std_last_%s' %value] = xc[-value:].std()\n        \n        X.loc[seg_id, 'mean_first_%s' %value] = xc[:value].mean()\n        X.loc[seg_id, 'mean_last_%s' %value] = xc[-value:].mean()\n       \n        X.loc[seg_id, 'min_first_%s' %value] = xc[:value].min()\n        X.loc[seg_id, 'min_last_%s' %value] = xc[-value:].min()\n        \n        X.loc[seg_id, 'max_first_%s' %value] = xc[:value].max()\n        X.loc[seg_id, 'max_last_%s' %value] = xc[-value:].max()\n        \n        #X.loc[seg_id, 'mean_change_rate_last_%s' %value] = np.mean(np.nonzero((np.diff(xc[-value:]) / xc[-value:][:-1]))[0])\n        #X.loc[seg_id, 'mean_change_rate_first_%s' %value] = np.mean(np.nonzero((np.diff(xc[:value]) / xc[:value][:-1]))[0])\n        \n    X.loc[seg_id, 'max_to_min'] = xc.max() / np.abs(xc.min())\n    X.loc[seg_id, 'max_min_diff'] = xc.max() - np.abs(xc.min())\n    X.loc[seg_id, 'sum'] = xc.sum()\n    X.loc[seg_id, 'mean_change_abs'] = np.mean(np.diff(xc))\n    \n    X.loc[seg_id, 'q95'] = np.quantile(xc, 0.95)\n    X.loc[seg_id, 'q80'] = np.quantile(xc, 0.80)\n    X.loc[seg_id, 'q05'] = np.quantile(xc, 0.05)\n    X.loc[seg_id, 'q20'] = np.quantile(xc, 0.20)\n    #interquartile rang\n    X.loc[seg_id, 'iqr'] = np.subtract(*np.percentile(xc, [75, 25]))\n    X.loc[seg_id, 'q999'] = np.quantile(xc,0.999)\n    X.loc[seg_id, 'q001'] = np.quantile(xc,0.001)\n    \n    X.loc[seg_id, 'trend'] = add_linear_trend(xc)\n  # recursive_sta_lta\n    df = 1/(np.sort(realFFT))\n    X.loc[seg_id, 'recursive_sta_lta1_mean'] = recursive_sta_lta(xc, 500, 10000).mean()\n    X.loc[seg_id, 'recursive_sta_lta2_mean'] = recursive_sta_lta(xc, 5000, 100000).mean()\n    #X.loc[seg_id, 'recursive_sta_lta3_mean'] = recursive_sta_lta(xc, df[0]*50, df[0]*100).mean()\n    #X.loc[seg_id, 'recursive_sta_lta4_mean'] = recursive_sta_lta(xc, df[1]*100, df[1]*200).mean()\n    \n    X.loc[seg_id, 'Hann_window_mean'] = (convolve(xc, hann(150), mode='same') / sum(hann(150))).mean()\n    X.loc[seg_id, 'Hilbert_mean'] = np.abs(hilbert(xc)).mean()\n    \n    # windowing\n    for winlen in [10, 100, 1000]:\n        x_rollwindow_std = xc.rolling(window=winlen, win_type='cosine').std().dropna().values\n        X.loc[seg_id, 'ave_roll_std_' + str(winlen)] = x_rollwindow_std.mean()    \n        X.loc[seg_id, 'std_roll_std_' + str(winlen)] = x_rollwindow_std.std()\n        X.loc[seg_id, 'max_roll_std_' + str(winlen)] = x_rollwindow_std.max()   \n        X.loc[seg_id, 'min_roll_std_' + str(winlen)] = x_rollwindow_std.min()    \n        X.loc[seg_id, 'q01_roll_std_' + str(winlen)] = np.quantile(x_rollwindow_std, 0.01)\n        X.loc[seg_id, 'q05_roll_std_' + str(winlen)] = np.quantile(x_rollwindow_std, 0.05)\n        X.loc[seg_id, 'q95_roll_std_' + str(winlen)] = np.quantile(x_rollwindow_std, 0.95)\n        X.loc[seg_id, 'q99_roll_std_' + str(winlen)] = np.quantile(x_rollwindow_std, 0.99)\n        X.loc[seg_id, 'ave_change_abs_roll_std_' + str(winlen)] = np.mean(np.diff(x_rollwindow_std))\n        \n        x_rollwindow_mean = xc.rolling(window=winlen, win_type='cosine').mean().dropna().values\n        X.loc[seg_id, 'ave_roll_mean_' + str(winlen)] = x_rollwindow_mean.mean()\n        X.loc[seg_id, 'std_roll_mean_' + str(winlen)] = x_rollwindow_mean.std()\n        X.loc[seg_id, 'max_roll_mean_' + str(winlen)] = x_rollwindow_mean.max()\n        X.loc[seg_id, 'min_roll_mean_' + str(winlen)] = x_rollwindow_mean.min()\n        X.loc[seg_id, 'q01_roll_mean_' + str(winlen)] = np.quantile(x_rollwindow_mean, 0.01)\n        X.loc[seg_id, 'q05_roll_mean_' + str(winlen)] = np.quantile(x_rollwindow_mean, 0.05)\n        X.loc[seg_id, 'q95_roll_mean_' + str(winlen)] = np.quantile(x_rollwindow_mean, 0.95)\n        X.loc[seg_id, 'q99_roll_mean_' + str(winlen)] = np.quantile(x_rollwindow_mean, 0.99)\n        X.loc[seg_id, 'ave_change_abs_roll_mean_' + str(winlen)] = np.mean(np.diff(x_rollwindow_mean))\n        \n    for winlen in [500, 1000]:\n        #winlen = randint(500, 10000)\n        X.loc[seg_id, 'Moving_average_%s_mean' %winlen] = xc.rolling(window=winlen).mean().mean(skipna=True)\n        X.loc[seg_id, 'MA_%s_std_mean' %winlen] = xc.rolling(window=winlen).std().mean()\n        X.loc[seg_id,'MA_%s_BB_high_mean' %winlen] = (X.loc[seg_id, 'Moving_average_%s_mean' %winlen] + \\\n                                                      2 * X.loc[seg_id, 'MA_%s_std_mean' %winlen]).mean()\n        X.loc[seg_id,'MA_%s_BB_low_mean' %winlen] = (X.loc[seg_id, 'Moving_average_%s_mean' %winlen] - \\\n                                                      2 * X.loc[seg_id, 'MA_%s_std_mean' %winlen]).mean()\n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T19:52:07.741683Z","iopub.execute_input":"2021-07-03T19:52:07.742052Z","iopub.status.idle":"2021-07-03T19:52:07.784816Z","shell.execute_reply.started":"2021-07-03T19:52:07.742022Z","shell.execute_reply":"2021-07-03T19:52:07.783503Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### 3.2. Process the Train Dataset\n","metadata":{}},{"cell_type":"code","source":"print('Train dataset has {} rows and {} columns.'.format(train_df.shape[0], train_df.shape[1]))\nprint(\"There are {} files in the test folder.\".format(len(os.listdir(os.path.join(INPATH, 'test' )))))\nprint('Each test segment has {} rows.'.format(\\\n    pd.read_csv(glob.glob(os.path.join(INPATH, 'test', '*.csv'))[0]).shape[0]))","metadata":{"execution":{"iopub.status.busy":"2021-07-03T19:52:07.787368Z","iopub.execute_input":"2021-07-03T19:52:07.787719Z","iopub.status.idle":"2021-07-03T19:52:08.239515Z","shell.execute_reply.started":"2021-07-03T19:52:07.787686Z","shell.execute_reply":"2021-07-03T19:52:08.238285Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_segments = int(np.round(train_df.shape[0] / 150000))\nprint(\"Number of segments: \", train_segments)\n\ntrain_X = pd.DataFrame(index=range(train_segments), dtype=np.float32)\ntrain_y = pd.DataFrame(index=range(train_segments), dtype=np.float32, columns=['time_to_failure'])","metadata":{"execution":{"iopub.status.busy":"2021-07-03T19:52:08.241285Z","iopub.execute_input":"2021-07-03T19:52:08.241594Z","iopub.status.idle":"2021-07-03T19:52:08.256085Z","shell.execute_reply.started":"2021-07-03T19:52:08.241565Z","shell.execute_reply":"2021-07-03T19:52:08.254693Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Feature engineering\n# Split the train data into segments of the same dimension as the test files.\nseg_length = 150000\n# Iterate over all segments\nfor seg_id in tqdm(range(train_segments)):\n    seg = train_df.iloc[seg_id*seg_length:(seg_id+1)*seg_length]\n    create_features(seg_id, seg, train_X)\n    train_y.loc[seg_id, 'time_to_failure'] = seg['time_to_failure'].values[-1]","metadata":{"execution":{"iopub.status.busy":"2021-07-03T19:52:08.258133Z","iopub.execute_input":"2021-07-03T19:52:08.258639Z","iopub.status.idle":"2021-07-03T20:47:57.634443Z","shell.execute_reply.started":"2021-07-03T19:52:08.258587Z","shell.execute_reply":"2021-07-03T20:47:57.633258Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Scale the train data\nscaler = StandardScaler()\nscaler.fit(train_X)\nscaled_train_X = pd.DataFrame(scaler.transform(train_X), columns=train_X.columns)\n#scaled_train_X.head(4)\nscaled_train_X.info()","metadata":{"execution":{"iopub.status.busy":"2021-07-03T20:47:57.636824Z","iopub.execute_input":"2021-07-03T20:47:57.637260Z","iopub.status.idle":"2021-07-03T20:47:57.695401Z","shell.execute_reply.started":"2021-07-03T20:47:57.637228Z","shell.execute_reply":"2021-07-03T20:47:57.694463Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"np.any(np.isnan(train_X))\nfor col in train_y:\n    if train_y[col].isnull().values.any():\n        print(col)\n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T20:47:57.697151Z","iopub.execute_input":"2021-07-03T20:47:57.697461Z","iopub.status.idle":"2021-07-03T20:47:57.706085Z","shell.execute_reply.started":"2021-07-03T20:47:57.697432Z","shell.execute_reply":"2021-07-03T20:47:57.704986Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### 3.3. Process the Test Dataset\n<a id='sec3.3'></a>","metadata":{}},{"cell_type":"code","source":"# Check the shape of the submission file\nsubmission = pd.read_csv(os.path.join(INPATH,'sample_submission.csv'))\nsubmission.shape","metadata":{"execution":{"iopub.status.busy":"2021-07-03T20:47:57.707377Z","iopub.execute_input":"2021-07-03T20:47:57.707683Z","iopub.status.idle":"2021-07-03T20:47:57.750207Z","shell.execute_reply.started":"2021-07-03T20:47:57.707654Z","shell.execute_reply":"2021-07-03T20:47:57.748963Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_test = pd.DataFrame(index=range(submission.shape[0]), dtype=np.float32)\n\n# Load files located in the test directory\ntest_files = glob.glob(os.path.join( INPATH, 'test', '*_*.*'))\n\n# Feature engineering for the test set\n# Iterate over all files in the test directory and\n# create the test dataframe with new features\nall_test_files = len(test_files)\n#all_test_files =600\nfor seg_id in tqdm(range(all_test_files)):\n    seg = pd.read_csv(test_files[seg_id])\n    create_features(seg_id, seg, X_test)\n    \n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T20:47:57.751908Z","iopub.execute_input":"2021-07-03T20:47:57.752388Z","iopub.status.idle":"2021-07-03T21:24:25.886927Z","shell.execute_reply.started":"2021-07-03T20:47:57.752339Z","shell.execute_reply":"2021-07-03T21:24:25.885598Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### 3.4 Data Selection","metadata":{}},{"cell_type":"code","source":"Y_train = train_y['time_to_failure'].values\n\nselection = list(scaled_train_X.columns)\nX_train = scaled_train_X[selection].values\n\nx_test = X_test[selection].values\n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:24:25.888984Z","iopub.execute_input":"2021-07-03T21:24:25.889414Z","iopub.status.idle":"2021-07-03T21:24:25.901108Z","shell.execute_reply.started":"2021-07-03T21:24:25.889366Z","shell.execute_reply":"2021-07-03T21:24:25.899797Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"scaled_train_X[selection].describe()","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:24:25.902815Z","iopub.execute_input":"2021-07-03T21:24:25.903247Z","iopub.status.idle":"2021-07-03T21:24:26.194398Z","shell.execute_reply.started":"2021-07-03T21:24:25.903201Z","shell.execute_reply":"2021-07-03T21:24:26.193469Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 4. Evaluating Some Models\n<a id='sec4'></a>","metadata":{}},{"cell_type":"code","source":"#Build the Models\nrandom_state = 1\n\n# Step 1: create a list containing all estimators with their default parameters\nmodel_list = [LinearRegression(), Ridge(), Lasso(),\n          KNeighborsRegressor(), DecisionTreeRegressor(),\n          RandomForestRegressor(), GradientBoostingRegressor(), \n          AdaBoostRegressor()]\n\n\n# Step 2: calculate the cross-validation mean and standard deviation for the estimators\ncv_mean, cv_std = [], []\n\n \nfor mdl in model_list: \n    print('-----------------------------', mdl)\n    cv = cross_val_score(mdl, X_train, y = Y_train, scoring='neg_mean_squared_error', cv = 7, n_jobs = -1)\n    \n    cv_mean.append(abs(cv.mean()))\n    cv_std.append(cv.std())\n\n        \n# Step 3: create a dataframe and plot  means with error bars\ncv_total = pd.DataFrame({'Algorithm': ['Linear Regression',  'Ridge', 'Lasso Regression',\n         'K Neighbors Regressor', 'Decision Tree Regressor', \n         'Random Forest Regressor', 'Gradient Boosting Regressor',\n         'Adaboost Regressor'],\n                         'CV-Means': cv_mean, \n                         'CV-Errors': cv_std})\n\nsns.barplot(x='CV-Means', y='Algorithm', data = cv_total, palette = 'Set1', orient = 'h',\\\n            **{'xerr': cv_std})\nplt.xlabel('Mean Squared Error')\nplt.title('Cross Validation Scores')\n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:24:26.196236Z","iopub.execute_input":"2021-07-03T21:24:26.196720Z","iopub.status.idle":"2021-07-03T21:26:27.789840Z","shell.execute_reply.started":"2021-07-03T21:24:26.196662Z","shell.execute_reply":"2021-07-03T21:26:27.788834Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"cv_total","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:26:27.791680Z","iopub.execute_input":"2021-07-03T21:26:27.792047Z","iopub.status.idle":"2021-07-03T21:26:27.806182Z","shell.execute_reply.started":"2021-07-03T21:26:27.792012Z","shell.execute_reply":"2021-07-03T21:26:27.805146Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 5. Hyperparameter Tuning\n<a id='sec5'></a>\nWe're optimising parameters for those models that have the highest accuracy.","metadata":{}},{"cell_type":"markdown","source":"#### Ridge","metadata":{}},{"cell_type":"code","source":"#Hyperparameter search for Ridge\nridge = Ridge()\nparam_grid = {'alpha': [0.05, 0.5, 0.7]\n             }\n\ncv = RepeatedKFold(n_splits=10, n_repeats=3, random_state=1)                                 \ngrid_ridge = GridSearchCV(ridge, param_grid = param_grid, cv=cv , scoring='neg_mean_squared_error', verbose = True, n_jobs = -1);\ngrid_ridge_result = grid_ridge.fit(X_train, Y_train);\n# summarize results for rf\nprint('-------------------------------------------\\n')\nprint(\"\\nBest Performance: MSE= %f using \\n%s\" % (grid_ridge_result.best_score_, grid_ridge_result.best_params_));\n\nmeans = grid_ridge_result.cv_results_['mean_test_score']\nstds = grid_ridge_result.cv_results_['std_test_score']\nparams = grid_ridge_result.cv_results_['params']\n#for mean, stdev, param in zip(means, stds, params):#\n#    print(\"%f (%f) with: %r\" % (mean, stdev, param))","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:26:27.807954Z","iopub.execute_input":"2021-07-03T21:26:27.808349Z","iopub.status.idle":"2021-07-03T21:26:28.376746Z","shell.execute_reply.started":"2021-07-03T21:26:27.808315Z","shell.execute_reply":"2021-07-03T21:26:28.375471Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### Random Forest Regressor","metadata":{}},{"cell_type":"code","source":"#Hyperparameter search for Randome Forest\n# https://scikit-learn.org/stable/modules/generated/sklearn.ensemble.RandomForestRegressor.html\nrf = RandomForestRegressor()\nparam_grid = {'n_estimators': [100, 200],\n              'bootstrap': [True],\n              'max_depth': [3, 5],\n              #'max_features': ['auto','sqrt'],\n              'min_samples_leaf': [2, 3],\n              'min_samples_split': [2, 3]}\ncv = RepeatedKFold(n_splits=10, n_repeats=3, random_state=random_state)                                   \ngrid_rf = GridSearchCV(rf, param_grid = param_grid, cv = 7, scoring='neg_mean_squared_error', verbose = True, n_jobs = -1);\ngrid_rf_result = grid_rf.fit(X_train, Y_train);\n# summarize results for rf\nprint('-------------------------------------------\\n')\nprint(\"\\nBest Performance : MSE= %f using \\n%s\" % (grid_rf_result.best_score_, grid_rf_result.best_params_));\n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:26:28.379030Z","iopub.execute_input":"2021-07-03T21:26:28.381821Z","iopub.status.idle":"2021-07-03T21:35:43.198956Z","shell.execute_reply.started":"2021-07-03T21:26:28.381750Z","shell.execute_reply":"2021-07-03T21:35:43.197589Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### AdaBoost Regressor","metadata":{}},{"cell_type":"code","source":"param_dist = {\n 'n_estimators': [50, 100],\n 'learning_rate' : [0.03,0.05,0.1,0.3],\n 'loss' : ['linear', 'square', 'exponential']\n }\n\n#cv = RepeatedKFold(n_splits=5, n_repeats=3, random_state=random_state)  \nrand_ada = RandomizedSearchCV(AdaBoostRegressor(), param_distributions = param_dist,\n                            cv=10, n_iter = 7, n_jobs=-1, scoring='neg_mean_squared_error')\nrand_ada_result = rand_ada.fit(X_train, Y_train)\n\nprint('-------------------------------------------\\n')\nprint(\"\\nBest Performance: MSE= %f using \\n%s\" % (rand_ada_result.best_score_, rand_ada_result.best_params_));\n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:35:43.200455Z","iopub.execute_input":"2021-07-03T21:35:43.200797Z","iopub.status.idle":"2021-07-03T21:38:44.950183Z","shell.execute_reply.started":"2021-07-03T21:35:43.200766Z","shell.execute_reply":"2021-07-03T21:38:44.949147Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 6. Building Final Models\n<a id='sec6'></a>","metadata":{}},{"cell_type":"code","source":"estimators = [('rf', grid_rf_result.best_estimator_),\\\n             ('Ridge', grid_ridge_result.best_estimator_),\n             ('AdaBoost', rand_ada_result.best_estimator_)]\n                \ntuned_voting = VotingRegressor(estimators = estimators, n_jobs = -1)\n\ntuned_voting.fit(X_train, Y_train)\n\n\ncv = RepeatedKFold(n_splits=5, n_repeats=3, random_state=random_state)  \n\ncv_voting = cross_val_score(tuned_voting, X_train, Y_train, cv = 10, scoring='neg_mean_squared_error')\n\nprint ('Tuned Models - Ensemble\\n-----------------------')\nprint ('Voting: {}%'.format(np.round(cv_voting.mean(), 2)))\n","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:44:25.401313Z","iopub.execute_input":"2021-07-03T21:44:25.401724Z","iopub.status.idle":"2021-07-03T21:46:29.520322Z","shell.execute_reply.started":"2021-07-03T21:44:25.401688Z","shell.execute_reply":"2021-07-03T21:46:29.518990Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 7. Prediction\n<a id='sec7'></a>\nI'll use the Ridge model for the final prediction. ","metadata":{}},{"cell_type":"code","source":"#y_pred_ensemble = tuned_voting.predict(x_test)\ny_pred_ridge = grid_ridge.predict(x_test)","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:54:47.893102Z","iopub.execute_input":"2021-07-03T21:54:47.893572Z","iopub.status.idle":"2021-07-03T21:54:47.905004Z","shell.execute_reply.started":"2021-07-03T21:54:47.893535Z","shell.execute_reply":"2021-07-03T21:54:47.903663Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 6.1. Submission","metadata":{}},{"cell_type":"code","source":"submission.time_to_failure = y_pred_ridge\nsubmission.to_csv('submission.csv',index=False)\nsubmission.head(4)","metadata":{"execution":{"iopub.status.busy":"2021-07-03T21:55:03.582739Z","iopub.execute_input":"2021-07-03T21:55:03.583297Z","iopub.status.idle":"2021-07-03T21:55:03.609990Z","shell.execute_reply.started":"2021-07-03T21:55:03.583263Z","shell.execute_reply":"2021-07-03T21:55:03.608793Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## References\n[Gabriel Preda](https://www.kaggle.com/gpreda/lanl-earthquake-eda-and-prediction)<br>\n[obspy](https://docs.obspy.org/_modules/obspy/signal/trigger.html)","metadata":{}}]}