{"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":"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\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","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**RNN with CNN feature extraction****","metadata":{}},{"cell_type":"markdown","source":"Forecasting 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.","metadata":{}},{"cell_type":"markdown","source":"We will address when the earthquake will take place. Specifically we’ll predict the time remaining before laboratory earthquakes occur from real-time seismic data.","metadata":{}},{"cell_type":"markdown","source":"Dataset Description\nthe dataset consists of seismic signals to predict the timing of laboratory earthquakes. The data comes from a well-known experimental set-up used to study earthquake physics. The acoustic_data input signal is used to predict the time remaining before the next laboratory earthquake (time_to_failure).\n\nThe training data is a single, continuous segment of experimental data. The test data consists of a folder containing many small segments. The data within each test file is continuous, but the test files do not represent a continuous segment of the experiment; thus, the predictions cannot be assumed to follow the same regular pattern seen in the training file.\n\nFor each seg_id in the test folder, you should predict a single time_to_failure corresponding to the time between the last row of the segment and the next laboratory earthquake.\n\nFile descriptions\n\ntrain.csv - A single, continuous training segment of experimental data.\n\ntest - A folder containing many small segments of test data.\n\nsample_sumbission.csv - A sample submission file in the correct format.","metadata":{}},{"cell_type":"markdown","source":"**Data fields**\n\nacoustic_data - the seismic signal\n\ntime_to_failure - the time (in seconds) until the next laboratory earthquake \n\nseg_id - the test segment ids for which predictions should be made (one prediction per segment)","metadata":{}},{"cell_type":"code","source":"pip install tensorflow","metadata":{"_kg_hide-input":true,"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%matplotlib inline\n\nfrom os import listdir, makedirs\nfrom os.path import isfile, join, basename, splitext, isfile, exists\n\nimport numpy as np\nimport pandas as pd\n\nfrom tqdm import tqdm_notebook\n\nimport tensorflow as tf\nimport keras.backend as K\n\nimport keras\nfrom keras.models import Sequential, Model\nfrom keras.layers import Dropout, Dense, Flatten, BatchNormalization\nfrom keras.layers import Convolution1D, ZeroPadding1D, MaxPooling1D, GlobalAveragePooling1D, GlobalMaxPooling1D\nfrom keras.layers import Concatenate, Average, Maximum, CuDNNLSTM, CuDNNGRU, Bidirectional, TimeDistributed\nfrom keras.callbacks import Callback, EarlyStopping, ModelCheckpoint\nfrom keras.engine.input_layer import Input\nfrom keras.models import load_model\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\npd.set_option('precision', 30)\nnp.set_printoptions(precision = 30)\n\nnp.random.seed(7723)\ntf.random.set_seed(1090)\n","metadata":{"execution":{"iopub.status.busy":"2022-10-22T02:57:54.74991Z","iopub.execute_input":"2022-10-22T02:57:54.750535Z","iopub.status.idle":"2022-10-22T02:58:00.351174Z","shell.execute_reply.started":"2022-10-22T02:57:54.7505Z","shell.execute_reply":"2022-10-22T02:58:00.350202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Memory saving function credit to https://www.kaggle.com/gemartin/load-data-reduce-memory-usage\ndef reduce_mem_usage(df):\n    \"\"\" iterate through all the columns of a dataframe and modify the data type\n        to reduce memory usage.        \n    \"\"\"\n    #start_mem = df.memory_usage().sum() / 1024**2\n    #print('Memory usage of dataframe is {:.2f} MB'.format(start_mem))\n\n    for col in df.columns:\n        col_type = df[col].dtype\n\n        if col_type != object:\n            c_min = df[col].min()\n            c_max = df[col].max()\n            if str(col_type)[:3] == 'int':\n                if c_min > np.iinfo(np.int8).min and c_max < np.iinfo(np.int8).max:\n                    df[col] = df[col].astype(np.int8)\n                elif c_min > np.iinfo(np.int16).min and c_max < np.iinfo(np.int16).max:\n                    df[col] = df[col].astype(np.int16)\n                elif c_min > np.iinfo(np.int32).min and c_max < np.iinfo(np.int32).max:\n                    df[col] = df[col].astype(np.int32)\n                elif c_min > np.iinfo(np.int64).min and c_max < np.iinfo(np.int64).max:\n                    df[col] = df[col].astype(np.int64)  \n            else:\n                if c_min > np.finfo(np.float16).min and c_max < np.finfo(np.float16).max:\n                    df[col] = df[col].astype(np.float16)\n                elif c_min > np.finfo(np.float32).min and c_max < np.finfo(np.float32).max:\n                    df[col] = df[col].astype(np.float32)\n                else:\n                    df[col] = df[col].astype(np.float64)\n\n    #end_mem = df.memory_usage().sum() / 1024**2\n    #print('Memory usage after optimization is: {:.2f} MB'.format(end_mem))\n    #print('Decreased by {:.1f}%'.format(100 * (start_mem - end_mem) / start_mem))\n\n    return df","metadata":{"execution":{"iopub.status.busy":"2022-10-22T02:58:03.976371Z","iopub.execute_input":"2022-10-22T02:58:03.977031Z","iopub.status.idle":"2022-10-22T02:58:03.99025Z","shell.execute_reply.started":"2022-10-22T02:58:03.976993Z","shell.execute_reply":"2022-10-22T02:58:03.989078Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\ntrain_df = pd.read_csv('../input/LANL-Earthquake-Prediction/train.csv', dtype={'acoustic_data': np.int8, 'time_to_failure': np.float32})","metadata":{"execution":{"iopub.status.busy":"2022-10-22T02:58:09.450163Z","iopub.execute_input":"2022-10-22T02:58:09.450528Z","iopub.status.idle":"2022-10-22T03:01:19.761173Z","shell.execute_reply.started":"2022-10-22T02:58:09.450496Z","shell.execute_reply":"2022-10-22T03:01:19.759967Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"I am using the earthquake data here. It is a time series data having two fields.\n\nacoustic_data = the time when signal was generated at the epicenter\ntime_to_failure = the time taken in seconds when quake is felt on the surface of earth","metadata":{}},{"cell_type":"code","source":"reduce_mem_usage(train_df)","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:01:19.763172Z","iopub.execute_input":"2022-10-22T03:01:19.763644Z","iopub.status.idle":"2022-10-22T03:01:27.873409Z","shell.execute_reply.started":"2022-10-22T03:01:19.763605Z","shell.execute_reply":"2022-10-22T03:01:27.872416Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_df.head()","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:01:59.646247Z","iopub.execute_input":"2022-10-22T03:01:59.646651Z","iopub.status.idle":"2022-10-22T03:01:59.657492Z","shell.execute_reply.started":"2022-10-22T03:01:59.646617Z","shell.execute_reply":"2022-10-22T03:01:59.656471Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train = train_df.acoustic_data.values\ny_train = train_df.time_to_failure.values","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:02:07.127223Z","iopub.execute_input":"2022-10-22T03:02:07.127626Z","iopub.status.idle":"2022-10-22T03:02:07.132377Z","shell.execute_reply.started":"2022-10-22T03:02:07.12759Z","shell.execute_reply":"2022-10-22T03:02:07.131204Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ends_mask = np.less(y_train[:-1], y_train[1:])\nsegment_ends = np.nonzero(ends_mask)\n\ntrain_segments = []\nstart = 0\nfor end in segment_ends[0]:\n    train_segments.append((start, end))\n    start = end\n    \nprint(train_segments)","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:02:10.037075Z","iopub.execute_input":"2022-10-22T03:02:10.037856Z","iopub.status.idle":"2022-10-22T03:02:15.422078Z","shell.execute_reply.started":"2022-10-22T03:02:10.0378Z","shell.execute_reply":"2022-10-22T03:02:15.421043Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.title('Segment sizes')\n_ = plt.bar(np.arange(len(train_segments)), [ s[1] - s[0] for s in train_segments])","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:02:19.076009Z","iopub.execute_input":"2022-10-22T03:02:19.076363Z","iopub.status.idle":"2022-10-22T03:02:19.339645Z","shell.execute_reply.started":"2022-10-22T03:02:19.076324Z","shell.execute_reply":"2022-10-22T03:02:19.338534Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class EarthQuakeRandom(keras.utils.all_utils.Sequence):\n\n    def __init__(self, x, y, x_mean, x_std, segments, ts_length, batch_size, steps_per_epoch):\n        self.x = x\n        self.y = y\n        self.segments = segments\n        self.ts_length = ts_length\n        self.batch_size = batch_size\n        self.steps_per_epoch = steps_per_epoch\n        self.segments_size = np.array([s[1] - s[0] for s in segments])\n        self.segments_p = self.segments_size / self.segments_size.sum()\n        self.x_mean = x_mean\n        self.x_std = x_std\n\n    def get_batch_size(self):\n        return self.batch_size\n\n    def get_ts_length(self):\n        return self.ts_length\n\n    def get_segments(self):\n        return self.segments\n\n    def get_segments_p(self):\n        return self.segments_p\n\n    def get_segments_size(self):\n        return self.segments_size\n\n    def __len__(self):\n        return self.steps_per_epoch\n\n    def __getitem__(self, idx):\n        segment_index = np.random.choice(range(len(self.segments)), p=self.segments_p)\n        segment = self.segments[segment_index]\n        end_indexes = np.random.randint(segment[0] + self.ts_length, segment[1], size=self.batch_size)\n\n        x_batch = np.empty((self.batch_size, self.ts_length))\n        y_batch = np.empty(self.batch_size, )\n\n        for i, end in enumerate(end_indexes):\n            x_batch[i, :] = self.x[end - self.ts_length: end]\n            y_batch[i] = self.y[end - 1]\n            \n        x_batch = (x_batch - self.x_mean)/self.x_std\n\n        return np.expand_dims(x_batch, axis=2), y_batch","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:02:22.59411Z","iopub.execute_input":"2022-10-22T03:02:22.594504Z","iopub.status.idle":"2022-10-22T03:02:22.608392Z","shell.execute_reply.started":"2022-10-22T03:02:22.594469Z","shell.execute_reply":"2022-10-22T03:02:22.607074Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"t_segments = [train_segments[i] for i in [ 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15]]\nv_segments = [train_segments[i] for i in [ 0, 1, 2, 3]]","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:02:27.93895Z","iopub.execute_input":"2022-10-22T03:02:27.939305Z","iopub.status.idle":"2022-10-22T03:02:27.944921Z","shell.execute_reply.started":"2022-10-22T03:02:27.939273Z","shell.execute_reply":"2022-10-22T03:02:27.943635Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Validating the data and calculate mean and standrad deviation on the training data.","metadata":{}},{"cell_type":"code","source":"x_sum = 0.\ncount = 0\n\nfor s in t_segments:\n    x_sum += X_train[s[0]:s[1]].sum()\n    count += (s[1] - s[0])\n\nX_train_mean = x_sum/count\n\nx2_sum = 0.\nfor s in t_segments:\n    x2_sum += np.power(X_train[s[0]:s[1]] - X_train_mean, 2).sum()\n\nX_train_std =  np.sqrt(x2_sum/count)\n\nprint(X_train_mean, X_train_std)\n","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:02:33.577402Z","iopub.execute_input":"2022-10-22T03:02:33.577782Z","iopub.status.idle":"2022-10-22T03:02:47.069964Z","shell.execute_reply.started":"2022-10-22T03:02:33.577749Z","shell.execute_reply":"2022-10-22T03:02:47.067982Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_gen = EarthQuakeRandom(\n    x = X_train, \n    y = y_train,\n    x_mean = X_train_mean, \n    x_std = X_train_std,\n    segments = t_segments,\n    ts_length = 150000,\n    batch_size = 64,\n    steps_per_epoch = 400\n)\n\nvalid_gen = EarthQuakeRandom(\n    x = X_train, \n    y = y_train,\n    x_mean = X_train_mean, \n    x_std = X_train_std,\n    segments = v_segments,\n    ts_length = 150000,\n    batch_size = 64,\n    steps_per_epoch = 400\n)","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:02:50.980894Z","iopub.execute_input":"2022-10-22T03:02:50.981843Z","iopub.status.idle":"2022-10-22T03:02:50.988588Z","shell.execute_reply.started":"2022-10-22T03:02:50.981789Z","shell.execute_reply":"2022-10-22T03:02:50.987261Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def CnnRnnModel():\n    i = Input(shape = (150000, 1))\n    \n    x = Convolution1D( 8, kernel_size = 10, strides = 10, activation='relu')(i)\n    x = Convolution1D(16, kernel_size = 10, strides = 10, activation='relu')(x)\n    x = Convolution1D(16, kernel_size = 10, strides = 10, activation='relu')(x)\n    x = CuDNNGRU(24, return_sequences = False, return_state = False)(x)\n    y = Dense(1)(x)\n\n    return Model(inputs = [i], outputs = [y])","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:02:56.007914Z","iopub.execute_input":"2022-10-22T03:02:56.008319Z","iopub.status.idle":"2022-10-22T03:02:56.015894Z","shell.execute_reply.started":"2022-10-22T03:02:56.008282Z","shell.execute_reply":"2022-10-22T03:02:56.014433Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model = CnnRnnModel()\nmodel.compile(loss='mean_absolute_error', optimizer='adam')\nmodel.summary()","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:03:00.370775Z","iopub.execute_input":"2022-10-22T03:03:00.371244Z","iopub.status.idle":"2022-10-22T03:03:03.273376Z","shell.execute_reply.started":"2022-10-22T03:03:00.371203Z","shell.execute_reply":"2022-10-22T03:03:03.272434Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"hist = model.fit_generator(\n    generator =  train_gen,\n    epochs = 50, \n    verbose = 0, \n    validation_data = valid_gen,\n    callbacks = [\n        EarlyStopping(monitor='val_loss', patience = 5, verbose = 1),\n        ModelCheckpoint(filepath='cnn_rnn.h5', monitor='val_loss', save_best_only=True, verbose=1)]\n)","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:03:12.553797Z","iopub.execute_input":"2022-10-22T03:03:12.554286Z","iopub.status.idle":"2022-10-22T03:22:55.812725Z","shell.execute_reply.started":"2022-10-22T03:03:12.554242Z","shell.execute_reply":"2022-10-22T03:22:55.811693Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.plot(hist.history['loss'])\nplt.plot(hist.history['val_loss'])\nplt.title('Model loss')\nplt.ylabel('Loss')\nplt.xlabel('Epoch')\n_= plt.legend(['Train', 'Test'], loc='upper left')","metadata":{"execution":{"iopub.status.busy":"2022-10-19T14:20:39.774682Z","iopub.execute_input":"2022-10-19T14:20:39.775079Z","iopub.status.idle":"2022-10-19T14:20:40.022483Z","shell.execute_reply.started":"2022-10-19T14:20:39.775047Z","shell.execute_reply":"2022-10-19T14:20:40.021416Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def load_test(ts_length = 150000):\n    base_dir = '../input/LANL-Earthquake-Prediction/test/'\n    test_files = [f for f in listdir(base_dir) if isfile(join(base_dir, f))]\n\n    ts = np.empty([len(test_files), ts_length])\n    ids = []\n    \n    i = 0\n    for f in tqdm_notebook(test_files):\n        ids.append(splitext(f)[0])\n        t_df = pd.read_csv(base_dir + f, dtype={\"acoustic_data\": np.int8})\n        ts[i, :] = t_df['acoustic_data'].values\n        i = i + 1\n\n    return ts, ids","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:24:08.73998Z","iopub.execute_input":"2022-10-22T03:24:08.741889Z","iopub.status.idle":"2022-10-22T03:24:08.753764Z","shell.execute_reply.started":"2022-10-22T03:24:08.741844Z","shell.execute_reply":"2022-10-22T03:24:08.752832Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_data, test_ids = load_test()","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:24:14.808158Z","iopub.execute_input":"2022-10-22T03:24:14.808891Z","iopub.status.idle":"2022-10-22T03:25:04.382626Z","shell.execute_reply.started":"2022-10-22T03:24:14.808848Z","shell.execute_reply":"2022-10-22T03:25:04.381623Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_test = ((test_data - X_train_mean)/ X_train_std).astype('float32')\nX_test = np.expand_dims(X_test, 2)\nX_test.shape","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:25:53.363199Z","iopub.execute_input":"2022-10-22T03:25:53.364073Z","iopub.status.idle":"2022-10-22T03:25:59.152495Z","shell.execute_reply.started":"2022-10-22T03:25:53.364037Z","shell.execute_reply":"2022-10-22T03:25:59.151442Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model = load_model('cnn_rnn.h5')","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:26:05.038156Z","iopub.execute_input":"2022-10-22T03:26:05.03914Z","iopub.status.idle":"2022-10-22T03:26:05.465125Z","shell.execute_reply.started":"2022-10-22T03:26:05.039094Z","shell.execute_reply":"2022-10-22T03:26:05.464123Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"y_pred = model.predict(X_test)","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:26:13.045664Z","iopub.execute_input":"2022-10-22T03:26:13.046504Z","iopub.status.idle":"2022-10-22T03:26:17.557601Z","shell.execute_reply.started":"2022-10-22T03:26:13.046439Z","shell.execute_reply":"2022-10-22T03:26:17.556601Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"names = ['x', 'y', 'z']\nindex = pd.MultiIndex.from_product([range(s)for s in X_test.shape], names=names)\ndf = pd.DataFrame({'A': X_test.flatten()}, index=index)['X_test']","metadata":{"execution":{"iopub.status.busy":"2022-10-22T03:28:32.266355Z","iopub.execute_input":"2022-10-22T03:28:32.266766Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!pip install xlsxwriter","metadata":{"execution":{"iopub.status.busy":"2022-10-19T14:25:35.397673Z","iopub.execute_input":"2022-10-19T14:25:35.398567Z","iopub.status.idle":"2022-10-19T14:25:46.870812Z","shell.execute_reply.started":"2022-10-19T14:25:35.398433Z","shell.execute_reply":"2022-10-19T14:25:46.869575Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"y_test_data = pd.DataFrame(X_test)\nwriter = pd.ExcelWriter(\"y-test-RNN.xlsx\", engine='xlsxwriter')\ny_test_data.to_excel(writer,sheet_name = \"sheet1_actual_data\", index=False)\nwriter.save()\n\ny_pred_data = pd.DataFrame(y_pred)\nwriter = pd.ExcelWriter(\"y_pred_RNN.xlsx\", engine='xlsxwriter')\ny_pred_data.to_excel(writer,sheet_name = \"sheet1_actual_data\", index=False)\nwriter.save()","metadata":{"execution":{"iopub.status.busy":"2022-10-19T14:25:46.945064Z","iopub.status.idle":"2022-10-19T14:25:46.945819Z","shell.execute_reply.started":"2022-10-19T14:25:46.945541Z","shell.execute_reply":"2022-10-19T14:25:46.945565Z"},"trusted":true},"execution_count":null,"outputs":[]}]}