{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport os\nimport psutil\n\nimport matplotlib.pyplot as plt\n%matplotlib inline\nfrom tqdm import tqdm_notebook\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.svm import NuSVR, SVR\nfrom sklearn.metrics import mean_absolute_error\npd.options.display.precision = 15\n\nimport lightgbm as lgb\nimport xgboost as xgb\nimport time\nimport datetime\n\nfrom sklearn.preprocessing import LabelEncoder\nfrom sklearn.model_selection import StratifiedKFold, KFold, RepeatedKFold\nfrom sklearn.metrics import mean_squared_error, mean_absolute_error\nfrom sklearn.linear_model import Ridge, RidgeCV\nimport gc\nfrom catboost import CatBoostRegressor\nimport seaborn as sns\nfrom scipy import stats","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":true},"cell_type":"code","source":"%%time\ntrain = pd.read_csv('../input/train.csv', dtype={'acoustic_data': np.int16, 'time_to_failure': np.float32})","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"2b369f9462919c0c6e29f8a5b6183cb3c57f057a","scrolled":true},"cell_type":"code","source":"rows = 150000\nsegments = int(np.floor(train.shape[0] / rows))\nX_tr = pd.DataFrame(index=range(segments), dtype=np.float64)\ny_tr = pd.DataFrame(index=range(segments), dtype=np.float64,\n                       columns=['time_to_failure'])\n\nfor segment in tqdm_notebook(range(segments)):\n    seg = train.iloc[segment*rows:segment*rows+rows]\n    x_raw = seg['acoustic_data']\n    x = x_raw.values\n    y = seg['time_to_failure'].values[-1]\n    \n    y_tr.loc[segment, 'time_to_failure'] = y\n    X_tr.loc[segment, 'ave'] = x.mean()\n    X_tr.loc[segment, 'std'] = x.std()\n    X_tr.loc[segment, 'max'] = x.max()\n    X_tr.loc[segment, 'min'] = x.min()\n    X_tr.loc[segment, 'q01'] = np.quantile(x,0.01)\n    X_tr.loc[segment, 'q05'] = np.quantile(x,0.05)\n    X_tr.loc[segment, 'q95'] = np.quantile(x,0.95)\n    X_tr.loc[segment, 'q99'] = np.quantile(x,0.99)\n    X_tr.loc[segment, 'abs_median'] = np.median(np.abs(x))\n    X_tr.loc[segment, 'abs_q95'] = np.quantile(np.abs(x),0.95)\n    X_tr.loc[segment, 'abs_q99'] = np.quantile(np.abs(x),0.99)\n    X_tr.loc[segment, 'F_test'], X_tr.loc[segment, 'p_test'] = stats.f_oneway(x[:30000],x[30000:60000],x[60000:90000],x[90000:120000],x[120000:])\n    X_tr.loc[segment, 'av_change_abs'] = np.mean(np.diff(x))\n    X_tr.loc[segment, 'av_change_rate'] = np.mean(np.nonzero((np.diff(x) / x[:-1]))[0])\n    X_tr.loc[segment, 'abs_max'] = np.abs(x).max()\n    \n    for windows in [10,100]:\n        x_roll_std = x_raw.rolling(windows).std().dropna().values\n        x_roll_mean = x_raw.rolling(windows).mean().dropna().values\n        \n        X_tr.loc[segment, 'ave_roll_std_' + str(windows)] = x_roll_std.mean()\n        X_tr.loc[segment, 'std_roll_std_' + str(windows)] = x_roll_std.std()\n        X_tr.loc[segment, 'max_roll_std_' + str(windows)] = x_roll_std.max()\n        X_tr.loc[segment, 'min_roll_std_' + str(windows)] = x_roll_std.min()\n        X_tr.loc[segment, 'q01_roll_std_' + str(windows)] = np.quantile(x_roll_std,0.01)\n        X_tr.loc[segment, 'q05_roll_std_' + str(windows)] = np.quantile(x_roll_std,0.05)\n        X_tr.loc[segment, 'q95_roll_std_' + str(windows)] = np.quantile(x_roll_std,0.95)\n        X_tr.loc[segment, 'q99_roll_std_' + str(windows)] = np.quantile(x_roll_std,0.99)\n        X_tr.loc[segment, 'av_change_abs_roll_std_' + str(windows)] = np.mean(np.diff(x_roll_std))\n        X_tr.loc[segment, 'av_change_rate_roll_std_' + str(windows)] = np.mean(np.nonzero((np.diff(x_roll_std) / x_roll_std[:-1]))[0])\n        X_tr.loc[segment, 'abs_max_roll_std_' + str(windows)] = np.abs(x_roll_std).max()\n        \n        X_tr.loc[segment, 'ave_roll_mean_' + str(windows)] = x_roll_mean.mean()\n        X_tr.loc[segment, 'std_roll_mean_' + str(windows)] = x_roll_mean.std()\n        X_tr.loc[segment, 'max_roll_mean_' + str(windows)] = x_roll_mean.max()\n        X_tr.loc[segment, 'min_roll_mean_' + str(windows)] = x_roll_mean.min()\n        X_tr.loc[segment, 'q01_roll_mean_' + str(windows)] = np.quantile(x_roll_mean,0.01)\n        X_tr.loc[segment, 'q05_roll_mean_' + str(windows)] = np.quantile(x_roll_mean,0.05)\n        X_tr.loc[segment, 'q95_roll_mean_' + str(windows)] = np.quantile(x_roll_mean,0.95)\n        X_tr.loc[segment, 'q99_roll_mean_' + str(windows)] = np.quantile(x_roll_mean,0.99)\n        X_tr.loc[segment, 'av_change_abs_roll_mean_' + str(windows)] = np.mean(np.diff(x_roll_mean))\n        X_tr.loc[segment, 'av_change_rate_roll_mean_' + str(windows)] = np.mean(np.nonzero((np.diff(x_roll_mean) / x_roll_mean[:-1]))[0])\n        X_tr.loc[segment, 'abs_max_roll_mean_' + str(windows)] = np.abs(x_roll_mean).max()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"287485e084611616118889d54754d8979b61449d"},"cell_type":"code","source":"X_tr.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"04367ea106ce4b11f53198981e9bde1d846265f2","scrolled":true},"cell_type":"code","source":"X_tr.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"8a8fefba8af693620aa3598ac3de8f2ad521bbf4"},"cell_type":"code","source":"submission = pd.read_csv('../input/sample_submission.csv', index_col='seg_id')\nX_test = pd.DataFrame(columns=X_tr.columns, dtype=np.float64, index=submission.index)\nplt.figure(figsize=(22, 16))\n\nfor i, seg_id in enumerate(tqdm_notebook(X_test.index)):\n    seg = pd.read_csv('../input/test/' + seg_id + '.csv')\n    \n    x_raw = seg['acoustic_data']\n    x_roll = x_raw.rolling(windows).std().dropna().values\n    x = x_raw.values\n    \n    X_test.loc[seg_id, 'ave'] = x.mean()\n    X_test.loc[seg_id, 'std'] = x.std()\n    X_test.loc[seg_id, 'max'] = x.max()\n    X_test.loc[seg_id, 'min'] = x.min()\n    X_test.loc[seg_id, 'q01'] = np.quantile(x,0.01)\n    X_test.loc[seg_id, 'q05'] = np.quantile(x,0.05)\n    X_test.loc[seg_id, 'q95'] = np.quantile(x,0.95)\n    X_test.loc[seg_id, 'q99'] = np.quantile(x,0.99)\n    X_test.loc[seg_id, 'abs_median'] = np.median(np.abs(x))\n    X_test.loc[seg_id, 'abs_q95'] = np.quantile(np.abs(x),0.95)\n    X_test.loc[seg_id, 'abs_q99'] = np.quantile(np.abs(x),0.99)\n    X_test.loc[seg_id, 'F_test'], X_tr.loc[segment, 'p_test'] = stats.f_oneway(x[:30000],x[30000:60000],x[60000:90000],x[90000:120000],x[120000:])\n    X_test.loc[seg_id, 'av_change_abs'] = np.mean(np.diff(x))\n    X_test.loc[seg_id, 'av_change_rate'] = np.mean(np.nonzero((np.diff(x) / x[:-1]))[0])\n    X_test.loc[seg_id, 'abs_max'] = np.abs(x).max()\n    \n    for windows in [10,100]:\n        x_roll_std = x_raw.rolling(windows).std().dropna().values\n        x_roll_mean = x_raw.rolling(windows).mean().dropna().values\n        \n        X_test.loc[seg_id, 'ave_roll_std_' + str(windows)] = x_roll_std.mean()\n        X_test.loc[seg_id, 'std_roll_std_' + str(windows)] = x_roll_std.std()\n        X_test.loc[seg_id, 'max_roll_std_' + str(windows)] = x_roll_std.max()\n        X_test.loc[seg_id, 'min_roll_std_' + str(windows)] = x_roll_std.min()\n        X_test.loc[seg_id, 'q01_roll_std_' + str(windows)] = np.quantile(x_roll_std,0.01)\n        X_test.loc[seg_id, 'q05_roll_std_' + str(windows)] = np.quantile(x_roll_std,0.05)\n        X_test.loc[seg_id, 'q95_roll_std_' + str(windows)] = np.quantile(x_roll_std,0.95)\n        X_test.loc[seg_id, 'q99_roll_std_' + str(windows)] = np.quantile(x_roll_std,0.99)\n        X_test.loc[seg_id, 'av_change_abs_roll_std_' + str(windows)] = np.mean(np.diff(x_roll_std))\n        X_test.loc[seg_id, 'av_change_rate_roll_std_' + str(windows)] = np.mean(np.nonzero((np.diff(x_roll_std) / x_roll_std[:-1]))[0])\n        X_test.loc[seg_id, 'abs_max_roll_std_' + str(windows)] = np.abs(x_roll_std).max()\n        \n        X_test.loc[seg_id, 'ave_roll_mean_' + str(windows)] = x_roll_mean.mean()\n        X_test.loc[seg_id, 'std_roll_mean_' + str(windows)] = x_roll_mean.std()\n        X_test.loc[seg_id, 'max_roll_mean_' + str(windows)] = x_roll_mean.max()\n        X_test.loc[seg_id, 'min_roll_mean_' + str(windows)] = x_roll_mean.min()\n        X_test.loc[seg_id, 'q01_roll_mean_' + str(windows)] = np.quantile(x_roll_mean,0.01)\n        X_test.loc[seg_id, 'q05_roll_mean_' + str(windows)] = np.quantile(x_roll_mean,0.05)\n        X_test.loc[seg_id, 'q95_roll_mean_' + str(windows)] = np.quantile(x_roll_mean,0.95)\n        X_test.loc[seg_id, 'q99_roll_mean_' + str(windows)] = np.quantile(x_roll_mean,0.99)\n        X_test.loc[seg_id, 'av_change_abs_roll_mean_' + str(windows)] = np.mean(np.diff(x_roll_mean))\n        X_test.loc[seg_id, 'av_change_rate_roll_mean_' + str(windows)] = np.mean(np.nonzero((np.diff(x_roll_mean) / x_roll_mean[:-1]))[0])\n        X_test.loc[seg_id, 'abs_max_roll_mean_' + str(windows)] = np.abs(x_roll_mean).max()\n        \n    \n    if i < 12:\n        plt.subplot(3, 4, i + 1)\n        plt.plot(seg['acoustic_data'])\n        plt.title(seg_id)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"9172a3903b178e151d8aefbc68e01426a39e9e2e"},"cell_type":"code","source":"n_fold = 5\nfolds = KFold(n_splits=n_fold, shuffle=True, random_state=11)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"815e13e6c0a78fc3fa461d4c823c97b185838a7e","_kg_hide-input":true},"cell_type":"code","source":"def train_model(X=X_tr, X_test=X_test, y=y_tr, params=None, folds=folds, model_type='lgb', plot_feature_importance=False, model=None):\n\n    oof = np.zeros(len(X))\n    prediction = np.zeros(len(X_test))\n    scores = []\n    feature_importance = pd.DataFrame()\n    for fold_n, (train_index, valid_index) in enumerate(folds.split(X)):\n        print('Fold', fold_n, 'started at', time.ctime())\n        X_train, X_valid = X.iloc[train_index], X.iloc[valid_index]\n        y_train, y_valid = y.iloc[train_index], y.iloc[valid_index]\n        \n        if model_type == 'lgb':\n            model = lgb.LGBMRegressor(**params, n_estimators = 20000, nthread = 4, n_jobs = -1)\n            model.fit(X_train, y_train, \n                    eval_set=[(X_train, y_train), (X_valid, y_valid)], eval_metric='mae',\n                    verbose=1000, early_stopping_rounds=200)\n            \n            y_pred_valid = model.predict(X_valid)\n            y_pred = model.predict(X_test, num_iteration=model.best_iteration_)\n            \n        if model_type == 'xgb':\n            train_data = xgb.DMatrix(data=X_train, label=y_train, feature_names=X_tr.columns)\n            valid_data = xgb.DMatrix(data=X_valid, label=y_valid, feature_names=X_tr.columns)\n\n            watchlist = [(train_data, 'train'), (valid_data, 'valid_data')]\n            model = xgb.train(dtrain=train_data, num_boost_round=20000, evals=watchlist, early_stopping_rounds=200, verbose_eval=500, params=params)\n            y_pred_valid = model.predict(xgb.DMatrix(X_valid, feature_names=X_tr.columns), ntree_limit=model.best_ntree_limit)\n            y_pred = model.predict(xgb.DMatrix(X_test, feature_names=X_tr.columns), ntree_limit=model.best_ntree_limit)\n            \n        if model_type == 'rcv':\n            model = RidgeCV(alphas=(0.01, 0.1, 1.0, 10.0, 100.0), scoring='neg_mean_absolute_error', cv=3)\n            model.fit(X_train, y_train)\n            print(model.alpha_)\n\n            y_pred_valid = model.predict(X_valid).reshape(-1,)\n            score = mean_absolute_error(y_valid, y_pred_valid)\n            print(f'Fold {fold_n}. MAE: {score:.4f}.')\n            print('')\n            \n            y_pred = model.predict(X_test).reshape(-1,)\n        \n        if model_type == 'sklearn':\n            model = model\n            model.fit(X_train, y_train)\n            \n            y_pred_valid = model.predict(X_valid).reshape(-1,)\n            score = mean_absolute_error(y_valid, y_pred_valid)\n            print(f'Fold {fold_n}. MAE: {score:.4f}.')\n            print('')\n            \n            y_pred = model.predict(X_test).reshape(-1,)\n        \n        if model_type == 'cat':\n            model = CatBoostRegressor(iterations=20000,  eval_metric='MAE', **params)\n            model.fit(X_train, y_train, eval_set=(X_valid, y_valid), cat_features=[], use_best_model=True, verbose=False)\n\n            y_pred_valid = model.predict(X_valid)\n            y_pred = model.predict(X_test)\n        \n        oof[valid_index] = y_pred_valid.reshape(-1,)\n        scores.append(mean_absolute_error(y_valid, y_pred_valid))\n\n        prediction += y_pred    \n        \n        if model_type == 'lgb':\n            # feature importance\n            fold_importance = pd.DataFrame()\n            fold_importance[\"feature\"] = X.columns\n            fold_importance[\"importance\"] = model.feature_importances_\n            fold_importance[\"fold\"] = fold_n + 1\n            feature_importance = pd.concat([feature_importance, fold_importance], axis=0)\n\n    prediction /= n_fold\n    \n    print('CV mean score: {0:.4f}.'.format(mean_absolute_error(y, oof)))\n    \n    if model_type == 'lgb':\n        feature_importance[\"importance\"] /= n_fold\n        if plot_feature_importance:\n            cols = feature_importance[[\"feature\", \"importance\"]].groupby(\"feature\").mean().sort_values(\n                by=\"importance\", ascending=False)[:50].index\n\n            best_features = feature_importance.loc[feature_importance.feature.isin(cols)]\n\n            plt.figure(figsize=(16, 12));\n            sns.barplot(x=\"importance\", y=\"feature\", data=best_features.sort_values(by=\"importance\", ascending=False));\n            plt.title('LGB Features (avg over folds)');\n        \n            return oof, prediction, feature_importance\n        return oof, prediction\n    \n    else:\n        return oof, prediction","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c532c79cb6f3417c3a9eec579979aaea91a0aa10"},"cell_type":"code","source":"params = {'num_leaves': 64,\n         'min_data_in_leaf': 50,\n         'objective': 'mae',\n         'max_depth': -1,\n         'learning_rate': 0.001,\n         \"boosting\": \"gbdt\",\n          \"feature_fraction\": 0.5,\n         \"bagging_freq\": 2,\n         \"bagging_fraction\": 0.5,\n         \"bagging_seed\": 0,\n         \"metric\": 'mae',\n         \"verbosity\": -1,\n         'reg_alpha': 1.0,\n         'reg_lambda': 1.0,\n         }\noof_lgb, prediction_lgb, feature_importance = train_model(params=params, model_type='lgb', plot_feature_importance=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"50d3d32bf0b7d71a1c4fe9adc00afeea06aa26f8"},"cell_type":"code","source":"params = {'num_leaves': 128,\n         'min_data_in_leaf': 50,\n         'objective': 'mae',\n         'max_depth': -1,\n         'learning_rate': 0.01,\n         \"boosting\": \"gbdt\",\n          \"feature_fraction\": 0.5,\n         \"bagging_freq\": 2,\n         \"bagging_fraction\": 0.7,\n         \"bagging_seed\": 0,\n         \"metric\": 'mae',\n         \"verbosity\": -1,\n         'reg_alpha': 2.0,\n         'reg_lambda': 1.5,\n         }\noof_lgb2, prediction_lgb2, feature_importance2 = train_model(params=params, model_type='lgb', plot_feature_importance=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1d35dd7aac6864e00d14b896a4ad551f23a864d6","scrolled":true},"cell_type":"code","source":"params = {'num_leaves': 32,\n         'min_data_in_leaf': 10,\n         'objective': 'mae',\n         'max_depth': -1,\n         'learning_rate': 0.001,\n         \"boosting\": \"gbdt\",\n          \"feature_fraction\": 0.8,\n         \"bagging_freq\": 2,\n         \"bagging_fraction\": 0.8,\n         \"bagging_seed\": 0,\n         \"metric\": 'mae',\n         \"verbosity\": -1,\n         'reg_alpha': 1.0,\n         'reg_lambda': 0.5,\n         }\noof_lgb3, prediction_lgb3, feature_importance3 = train_model(params=params, model_type='lgb', plot_feature_importance=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"2d01a97d68c4031b3b2d11974a17b2e7053efc88","scrolled":false},"cell_type":"code","source":"xgb_params = {'eta': 0.01, 'max_depth': 6, 'subsample': 0.8, 'colsample_bytree': 0.8,'colsample_bylevel': 0.8, 'colsample_bynode': 0.8,\n          'lambda': 0.1, 'alpha' : 0.1, 'objective': 'reg:linear', 'eval_metric': 'mae', 'silent': True, 'nthread': 4}\noof_xgb, prediction_xgb = train_model(params=xgb_params, model_type='xgb')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"695305ee9716cc7bdbbd3866388d24df01f58432","scrolled":true},"cell_type":"code","source":"xgb_params = {'eta': 0.01, 'max_depth': 12, 'subsample': 0.5, 'colsample_bytree': 0.5,'colsample_bylevel': 0.5, 'colsample_bynode': 0.5,\n              'lambda': 0.5, 'alpha' : 0.5, 'objective': 'reg:linear', 'eval_metric': 'mae', 'silent': True, 'nthread': 4}\noof_xgb2, prediction_xgb2 = train_model(params=xgb_params, model_type='xgb')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"6ed44e58c888d3891ece9bb8d01a5d56ae682eb1","scrolled":false},"cell_type":"code","source":"print(mean_absolute_error(y_tr, (oof_lgb + oof_lgb2  + oof_lgb3 + oof_xgb + oof_xgb2) /5))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"8359a6fd7cc8aec99fbad3b83636da4a47188e85"},"cell_type":"code","source":"plt.figure(figsize=(12, 8))\nplt.plot(y_tr, color='g', label='y_train')\nplt.plot(oof_lgb, color='b', label='lgb')\nplt.plot(oof_lgb2, color='b', label='lgb2')\nplt.plot(oof_lgb3, color='b', label='lgb3')\nplt.plot(oof_xgb, color='teal', label='xgb')\nplt.plot(oof_xgb2, color='teal', label='xgb')\nplt.plot((oof_lgb + oof_lgb2  + oof_lgb3 + oof_xgb + oof_xgb2) /5 , color='gold', label='blend')\nplt.legend();\nplt.title('Predictions vs actual');","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"aed2c67cf5b7596b785bb85b0943ee9fef6e8187"},"cell_type":"code","source":"prediction_lgb[:10], prediction_lgb2[:10], prediction_lgb3[:10], prediction_xgb[:10], prediction_xgb2[:10]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"acf417f3540ca69267b1088e199e4ae358b36439"},"cell_type":"code","source":"submission['time_to_failure'] = (prediction_lgb + prediction_lgb2 + prediction_lgb3 + prediction_xgb + prediction_xgb2) / 5\nprint(submission.head())\nsubmission.to_csv('submission_all.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"8c352806a79c1becb977292d10142d21a41db69a"},"cell_type":"code","source":"","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}