{"cells":[{"metadata":{"_uuid":"c2e7ee1afd50b46c9840880c2fabc8088374c41e"},"cell_type":"markdown","source":"# General information\nForecasting earthquakes is one of the most important problems in Earth science because of their devastating consequences. Current scientific studies related to earthquake forecasting focus on three key points: **when** the event will occur, **where** it will occur, and **how** large it will be. In this competition we try to predict time left to the next laboratory earthquake based on seismic signal data to answer the question of **when** earthquake will occur.\n\nTraining data represents one huge signal, but in test data we have many separate chunks, for each of which we need to predict time to failure."},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"scrolled":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 in \n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport seaborn as sns\nimport matplotlib.pyplot as plt\n\n%matplotlib inline\n\n# Input data files are available in the \"../input/\" directory.\n# For example, running this (by clicking run or pressing Shift+Enter) will list the files in the input directory\n\nimport os\nprint(os.listdir(\"../input\"))\n\n# Any results you write to the current directory are saved as output.","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","collapsed":true,"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":false},"cell_type":"markdown","source":"The data is huge, training data contains nearly **600 million rows** and that is A LOT of data to understand. I am just curious to know how long Kaggle servers will take to read this data."},{"metadata":{"trusted":true,"_uuid":"7bddcb9e3eceb1fad0ac8a4e206a0a6e8cc6a2ad"},"cell_type":"code","source":"%%time\ntrain = pd.read_csv('../input/train.csv', dtype={'acoustic_data': np.int16, 'time_to_failure': np.float64})","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"814c26a06ca23f8565342928bc378aac9d8b4df6"},"cell_type":"code","source":"train.head(5)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"afd59ef1a120142a44fe88f80c4ed81ad8d4daa2"},"cell_type":"code","source":"train.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"8395b5f344f51b69b746f41c0093853b826c2902"},"cell_type":"code","source":"partial_train = train[::20]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7020ae987fc8de156bcfec2d988a37090f996070"},"cell_type":"code","source":"partial_train.head(5)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"797739323dafa5bfee8b8e263f0983de25eb0984"},"cell_type":"code","source":"partial_train['time_to_failure'].shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1695f7c95dc35e12bd18cc247ec3009a5a6018fd"},"cell_type":"code","source":"figure, axes1 = plt.subplots(figsize=(18,10))\n\nplt.title(\"Seismic Data Trends with 5% sample of original data\")\n\nplt.plot(partial_train['acoustic_data'], color='r')\naxes1.set_ylabel('Acoustic Data', color='r')\nplt.legend(['Acoustic Data'])\n\naxes2 = axes1.twinx()\nplt.plot(partial_train['time_to_failure'], color='g')\naxes2.set_ylabel('Time to Failure', color='g')\nplt.legend(['Time to Failure'])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"f4cfd65d0f81ccce96b353248c7ec88b7317c6d6"},"cell_type":"markdown","source":"## Summary\nData is **HUGE!** <br>\nBefore every Failure, **there is a peak** in Acoustic Data. <br>\nSince we only have acoustic data to predict the time to failure, we need to **generate some features**. <br>\nI would use **Feature Engineering** to generate some common statistical features like **Mean, Variance, Max, Min, Std. Dev.** of our acoustic data. <br>\nSince the test data is **segmented into chunks of data**, it is better to segment our training data into chunks and then generate the features. <br>\n"},{"metadata":{"_uuid":"958887ad02712c6b4f2a5087bd4342deab475ce5"},"cell_type":"markdown","source":"# Data Preparation\nIt took about 6 million datapoints to get one failure. Each of test csv files contains 150,000 datapoints. Thus, we should be able to get about **629145480 / 150000 = 4194** samples (subsets similar to our test dataset for training our model) before the first failure in our training set. Their corresponding labels should be the time_to_failure at the last datapoint in each subset. We would use a regression model to predict the time to failure for each corresponding subset sample."},{"metadata":{"trusted":true,"_uuid":"b278aaed0dcf6f3f02e24a0e421a58f09044e427"},"cell_type":"code","source":"# list of features to be engineered\n\nfeatures = ['mean','max','variance','min', 'stdev', 'max-min-diff',\n            'max-mean-diff', 'mean-change-abs', 'abs-max', 'abs-min',\n            'std-first-50000', 'std-last-50000', 'mean-first-50000',\n            'mean-last-50000', 'max-first-50000', 'max-last-50000',\n            'min-first-50000', 'min-last-50000']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"498f357d68063891738fe0befbb66a00b782a7e7"},"cell_type":"code","source":"# Feature Engineering\n\nrows = 150000\nsegments = int(np.floor(train.shape[0] / rows))\n\nX = pd.DataFrame(index=range(segments), dtype=np.float64,\n                       columns=features)\nY = pd.DataFrame(index=range(segments), dtype=np.float64,\n                       columns=['time_to_failure'])\n\nfor segment in range(segments):\n    seg = train.iloc[segment*rows:segment*rows+rows]\n    x = seg['acoustic_data'].values\n    y = seg['time_to_failure'].values[-1]\n    \n    Y.loc[segment, 'time_to_failure'] = y\n    \n    X.loc[segment, 'mean'] = x.mean()\n    X.loc[segment, 'stdev'] = x.std()\n    X.loc[segment, 'variance'] = np.var(x)\n    X.loc[segment, 'max'] = x.max()\n    X.loc[segment, 'min'] = x.min()\n    X.loc[segment, 'max-min-diff'] = x.max()-x.min()\n    X.loc[segment, 'max-mean-diff'] = x.max()-x.mean()\n    X.loc[segment, 'mean-change-abs'] = np.mean(np.diff(x))\n    X.loc[segment, 'abs-min'] = np.abs(x).min()\n    X.loc[segment, 'abs-max'] = np.abs(x).max()\n    X.loc[segment, 'std-first-50000'] = x[:50000].std()\n    X.loc[segment, 'std-last-50000'] = x[-50000:].std()\n    X.loc[segment, 'mean-first-50000'] = x[:50000].min()\n    X.loc[segment, 'mean-last-50000'] = x[-50000:].mean()\n    X.loc[segment, 'max-first-50000'] = x[:50000].max()\n    X.loc[segment, 'max-last-50000'] = x[-50000:].max()\n    X.loc[segment, 'min-first-50000'] = x[:50000].min()\n    X.loc[segment, 'min-last-50000'] = x[-50000:].min()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a953de94af17fb11fe77bbf7bd980c16cfb9a352"},"cell_type":"code","source":"X.head(5)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4e699c228a3edffe4c3b5ac7a709a4b0834eeac7"},"cell_type":"code","source":"data = pd.concat([X,Y],axis=1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c64cd1871bf7c69b5e2a59326bdfbf34af197ee3"},"cell_type":"code","source":"sns.set(rc={'figure.figsize': (18,12)})\nsns.pairplot(data)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"2503064301968547ed69e76acace7e8500a6594b"},"cell_type":"code","source":"from sklearn.model_selection import train_test_split\nfrom sklearn.metrics import r2_score\nfrom sklearn.linear_model import LinearRegression\nfrom sklearn.metrics import mean_squared_error, mean_absolute_error\nfrom sklearn.preprocessing import StandardScaler","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"2df8175c20c7a3e54c80a2a39e03944bdd2b4352"},"cell_type":"code","source":"X_train,X_test, y_train, y_test = train_test_split(X, Y, test_size = 0.2, random_state=1210)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"227edacf38279c210f309403322c80d380ea6c5d"},"cell_type":"code","source":"scaler = StandardScaler()\nscaler.fit(X_train)\nX_train_sc = pd.DataFrame(scaler.transform(X_train), columns=X_train.columns)\n\nscaler.fit(X_test)\nX_test_sc = pd.DataFrame(scaler.transform(X_test), columns=X_test.columns)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e4f608cf79917939c271d2676c02423c6f5ecfaf"},"cell_type":"code","source":"lr = LinearRegression()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"b67285fdba042a3e04a9d1c8f2a1b81811675035"},"cell_type":"code","source":"lr.fit(X_train_sc,y_train)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c7edb73b779b8edbea7564851f3e6398c772ecd3"},"cell_type":"code","source":"pred = lr.predict(X_test_sc)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4116f8a545d505391994215226d14ce457e2e4a6"},"cell_type":"code","source":"mean_absolute_error(y_test, pred)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"2b08fd145663a866022f9596dd96003427910f06"},"cell_type":"code","source":"from lightgbm import LGBMRegressor\nparams = {'num_leaves': 54,'min_data_in_leaf': 79,'objective': 'huber',\n         'max_depth': -1, 'learning_rate': 0.01, \"boosting\": \"gbdt\",\n         # \"feature_fraction\": 0.8354507676881442,\n         \"bagging_freq\": 3,\"bagging_fraction\": 0.8126672064208567,\n         \"bagging_seed\": 11,\"metric\": 'mae',\n         \"verbosity\": -1,'reg_alpha': 1.1302650970728192,\n         'reg_lambda': 0.3603427518866501}","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"ab4b523f04d06804446ccae0b93ab6b8e419ac77"},"cell_type":"code","source":"lgbm = LGBMRegressor(nthread=4,n_estimators=10000,\n            learning_rate=0.01,num_leaves=54,\n            colsample_bytree=0.9497036,subsample=0.8715623,\n            max_depth=8,reg_alpha=0.04,\n            reg_lambda=0.073,min_child_weight=40,silent=-1,verbose=-1,)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"fbbe27ceb34aaf994190345eb5efe514dc9fc086"},"cell_type":"code","source":"lgbm.fit(X_train, y_train)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"889ced6e1dfc2197e6cfaadfd98e2ac1b1857800"},"cell_type":"code","source":"pred_lgbm = lgbm.predict(X_test)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"0bb169bce606638a8786497fb874b64b6917baca"},"cell_type":"code","source":"mean_absolute_error(y_test,pred_lgbm)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"3b394c9b9092d987623d43e01caec02e9abfa401"},"cell_type":"markdown","source":"The score can definitely be improved a lot. I will come back with more insights and update this kernel for better scores. Maybe using multiple predictors and the averaging their results will help reduce the error."}],"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}