{"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":"!pip3 install basemap","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:14.076846Z","iopub.execute_input":"2022-04-06T03:40:14.077509Z","iopub.status.idle":"2022-04-06T03:40:42.587705Z","shell.execute_reply.started":"2022-04-06T03:40:14.077418Z","shell.execute_reply":"2022-04-06T03:40:42.586840Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np \nimport pandas as pd\nimport matplotlib.pyplot as plt\nfrom mpl_toolkits.basemap import Basemap\nimport seaborn as sns\nsns.set(style=\"darkgrid\")","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:42.589544Z","iopub.execute_input":"2022-04-06T03:40:42.589844Z","iopub.status.idle":"2022-04-06T03:40:43.662825Z","shell.execute_reply.started":"2022-04-06T03:40:42.589800Z","shell.execute_reply":"2022-04-06T03:40:43.662136Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data = pd.read_csv('../input/earthquake-database/database.csv')","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:43.664218Z","iopub.execute_input":"2022-04-06T03:40:43.664488Z","iopub.status.idle":"2022-04-06T03:40:43.762903Z","shell.execute_reply.started":"2022-04-06T03:40:43.664451Z","shell.execute_reply":"2022-04-06T03:40:43.762173Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data.head()","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:43.764992Z","iopub.execute_input":"2022-04-06T03:40:43.765263Z","iopub.status.idle":"2022-04-06T03:40:43.799550Z","shell.execute_reply.started":"2022-04-06T03:40:43.765228Z","shell.execute_reply":"2022-04-06T03:40:43.798766Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data.shape","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:43.801149Z","iopub.execute_input":"2022-04-06T03:40:43.801425Z","iopub.status.idle":"2022-04-06T03:40:43.806661Z","shell.execute_reply.started":"2022-04-06T03:40:43.801388Z","shell.execute_reply":"2022-04-06T03:40:43.805839Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(\"Min Value: \"+ str(data['Magnitude'].min()))\nprint(\"Max Value: \" + str(data['Magnitude'].max()))","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:43.808874Z","iopub.execute_input":"2022-04-06T03:40:43.809209Z","iopub.status.idle":"2022-04-06T03:40:43.819581Z","shell.execute_reply.started":"2022-04-06T03:40:43.809172Z","shell.execute_reply":"2022-04-06T03:40:43.818848Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"g8 = data[data['Magnitude'] > 8]\ng8['Location Source'].value_counts()","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:43.820907Z","iopub.execute_input":"2022-04-06T03:40:43.821446Z","iopub.status.idle":"2022-04-06T03:40:43.832588Z","shell.execute_reply.started":"2022-04-06T03:40:43.821376Z","shell.execute_reply":"2022-04-06T03:40:43.831798Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.hist(data['Magnitude'])\n\nplt.xlabel('Magnitude Size')\nplt.ylabel('Number of Occurrences')","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:43.833811Z","iopub.execute_input":"2022-04-06T03:40:43.834113Z","iopub.status.idle":"2022-04-06T03:40:44.126831Z","shell.execute_reply.started":"2022-04-06T03:40:43.834076Z","shell.execute_reply":"2022-04-06T03:40:44.126163Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.countplot(x=\"Magnitude Type\", data=data)\nplt.ylabel('Frequency')\nplt.title('Magnitude Type VS Frequency')\nprint(\" local magnitude (ML), surface-wave magnitude (Ms), body-wave magnitude (Mb), moment magnitude (Mw)\")","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:44.127895Z","iopub.execute_input":"2022-04-06T03:40:44.128462Z","iopub.status.idle":"2022-04-06T03:40:44.416054Z","shell.execute_reply.started":"2022-04-06T03:40:44.128422Z","shell.execute_reply":"2022-04-06T03:40:44.415403Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def get_marker_color(magnitude):\n    if magnitude < 6.2:\n        return ('go')\n    elif magnitude < 7.5:\n        return ('yo')\n    else:\n        return ('ro')\n\nplt.figure(figsize=(14,10))\n\neq_map = Basemap(projection='robin', resolution = 'l',\n              lat_0=0, lon_0=-130)\neq_map.drawcoastlines()\neq_map.drawcountries()\neq_map.fillcontinents(color = 'gray')\neq_map.drawmapboundary()\neq_map.drawmeridians(np.arange(0, 360, 30))\neq_map.drawparallels(np.arange(-90, 90, 30))\n \n# read longitude, latitude and magnitude\nlons = data['Longitude'].values\nlats = data['Latitude'].values\nmagnitudes = data['Magnitude'].values\ntimestrings = data['Date'].tolist()\n    \nmin_marker_size = 0.5\nfor lon, lat, mag in zip(lons, lats, magnitudes):\n    x,y = eq_map(lon, lat)\n    msize = mag # * min_marker_size\n    marker_string = get_marker_color(mag)\n    eq_map.plot(x, y, marker_string, markersize=msize)\n    \ntitle_string = \"Earthquakes of Magnitude 5.5 or Greater\\n\"\ntitle_string += \"%s - %s\" % (timestrings[0][:10], timestrings[-1][:10])\nplt.title(title_string)\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:40:44.419515Z","iopub.execute_input":"2022-04-06T03:40:44.419724Z","iopub.status.idle":"2022-04-06T03:42:07.616791Z","shell.execute_reply.started":"2022-04-06T03:40:44.419698Z","shell.execute_reply":"2022-04-06T03:42:07.614910Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import datetime\ndata['date'] = data['Date'].apply(lambda x: pd.to_datetime(x))\ndata['year'] = data['date'].apply(lambda x: str(x).split('-')[0])\nplt.figure(figsize=(15, 8))\nsns.set(font_scale=1.0)\nsns.countplot(x=\"year\", data=data)\nplt.ylabel('Number Of Earthquakes')\nplt.title('Number of Earthquakes In Each Year')","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:42:07.618146Z","iopub.execute_input":"2022-04-06T03:42:07.618441Z","iopub.status.idle":"2022-04-06T03:42:11.322083Z","shell.execute_reply.started":"2022-04-06T03:42:07.618399Z","shell.execute_reply":"2022-04-06T03:42:11.321422Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data['year'].value_counts()[:1]","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:42:11.323143Z","iopub.execute_input":"2022-04-06T03:42:11.324509Z","iopub.status.idle":"2022-04-06T03:42:11.337517Z","shell.execute_reply.started":"2022-04-06T03:42:11.324468Z","shell.execute_reply":"2022-04-06T03:42:11.336731Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"x = data['year'].unique()\ny = data['year'].value_counts()\n\ncount = []\nfor i in range(len(x)):\n    key = x[i]\n    count.append(y[key])\n\nplt.figure(figsize=(10, 8))\n\nplt.scatter(x, count)\nplt.xlabel('Year')\nplt.ylabel('Number of Earthquakes')\nplt.title('Earthquakes Per year from 1995 to 2016')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:42:11.338896Z","iopub.execute_input":"2022-04-06T03:42:11.339349Z","iopub.status.idle":"2022-04-06T03:42:12.406544Z","shell.execute_reply.started":"2022-04-06T03:42:11.339313Z","shell.execute_reply":"2022-04-06T03:42:12.405662Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Magnitude Classes\n\n- **Disastrous**:   M > =8\n- **Major**:   7 < =M < 7.9\n- **Strong**:  6 < = M < 6.9\n- **Moderate**: 5.5 < =M < 5.9","metadata":{}},{"cell_type":"code","source":"data.loc[data['Magnitude'] >=8, 'Class'] = 'Disastrous'\ndata.loc[ (data['Magnitude'] >= 7) & (data['Magnitude'] < 7.9), 'Class'] = 'Major'\ndata.loc[ (data['Magnitude'] >= 6) & (data['Magnitude'] < 6.9), 'Class'] = 'Strong'\ndata.loc[ (data['Magnitude'] >= 5.5) & (data['Magnitude'] < 5.9), 'Class'] = 'Moderate'","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:42:12.407937Z","iopub.execute_input":"2022-04-06T03:42:12.408346Z","iopub.status.idle":"2022-04-06T03:42:12.429081Z","shell.execute_reply.started":"2022-04-06T03:42:12.408291Z","shell.execute_reply":"2022-04-06T03:42:12.428234Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Magnitude Class distribution\n\nsns.countplot(x=\"Class\", data=data)\nplt.ylabel('Frequency')\nplt.title('Magnitude Class VS Frequency')","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:42:12.430927Z","iopub.execute_input":"2022-04-06T03:42:12.431617Z","iopub.status.idle":"2022-04-06T03:42:12.776554Z","shell.execute_reply.started":"2022-04-06T03:42:12.431574Z","shell.execute_reply":"2022-04-06T03:42:12.775848Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np \nimport pandas as pd\nimport os\nfrom tqdm import tqdm","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:42:12.780229Z","iopub.execute_input":"2022-04-06T03:42:12.782181Z","iopub.status.idle":"2022-04-06T03:42:12.787610Z","shell.execute_reply.started":"2022-04-06T03:42:12.782142Z","shell.execute_reply":"2022-04-06T03:42:12.786920Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Fix seeds\nfrom numpy.random import seed\nseed(639)\nfrom tensorflow.random import set_seed\nset_seed(5944)","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:42:12.791683Z","iopub.execute_input":"2022-04-06T03:42:12.793963Z","iopub.status.idle":"2022-04-06T03:42:16.877543Z","shell.execute_reply.started":"2022-04-06T03:42:12.793920Z","shell.execute_reply":"2022-04-06T03:42:16.876736Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Import\nfloat_data = pd.read_csv(\"../input/LANL-Earthquake-Prediction/train.csv\", dtype={\"acoustic_data\": np.float32, \"time_to_failure\": np.float32}).values\n","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:42:16.879000Z","iopub.execute_input":"2022-04-06T03:42:16.879246Z","iopub.status.idle":"2022-04-06T03:45:38.200518Z","shell.execute_reply.started":"2022-04-06T03:42:16.879210Z","shell.execute_reply":"2022-04-06T03:45:38.199735Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Helper function for the data generator. Extracts mean, standard deviation, and quantiles per time step.\n# Can easily be extended. Expects a two dimensional array.\ndef extract_features(z):\n     return np.c_[z.mean(axis=1), \n                  z.min(axis=1),\n                  z.max(axis=1),\n                  z.std(axis=1)]","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:38.201950Z","iopub.execute_input":"2022-04-06T03:45:38.202199Z","iopub.status.idle":"2022-04-06T03:45:38.206736Z","shell.execute_reply.started":"2022-04-06T03:45:38.202164Z","shell.execute_reply":"2022-04-06T03:45:38.206094Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def create_X(x, last_index=None, n_steps=150, step_length=1000):\n    if last_index == None:\n        last_index=len(x)\n       \n    assert last_index - n_steps * step_length >= 0\n\n    # Reshaping and approximate standardization with mean 5 and std 3.\n    temp = (x[(last_index - n_steps * step_length):last_index].reshape(n_steps, -1) - 5 ) / 3\n    \n    # Extracts features of sequences of full length 1000, of the last 100 values and finally also \n    # of the last 10 observations. \n    return np.c_[extract_features(temp),\n                 extract_features(temp[:, -step_length // 10:]),\n                 extract_features(temp[:, -step_length // 100:])]","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:38.208297Z","iopub.execute_input":"2022-04-06T03:45:38.208831Z","iopub.status.idle":"2022-04-06T03:45:38.218786Z","shell.execute_reply.started":"2022-04-06T03:45:38.208795Z","shell.execute_reply":"2022-04-06T03:45:38.217895Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Query \"create_X\" to figure out the number of features\nn_features = create_X(float_data[0:150000]).shape[1]\nprint(\"Our RNN is based on %i features\"% n_features)","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:38.221699Z","iopub.execute_input":"2022-04-06T03:45:38.222010Z","iopub.status.idle":"2022-04-06T03:45:38.238582Z","shell.execute_reply.started":"2022-04-06T03:45:38.221983Z","shell.execute_reply":"2022-04-06T03:45:38.237841Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# The generator endlessly selects \"batch_size\" ending positions of sub-time series. For each ending position,\n# the \"time_to_failure\" serves as target, while the features are created by the function \"create_X\".\n\ndef generator(data, min_index=0, max_index=None, batch_size=16, n_steps=150, step_length=1000):\n    if max_index is None:\n        max_index = len(data) - 1\n     \n    while True:\n        # Pick indices of ending positions\n        rows = np.random.randint(min_index + n_steps * step_length, max_index, size=batch_size)\n         \n        # Initialize feature matrices and targets\n        samples = np.zeros((batch_size, n_steps, n_features))\n        targets = np.zeros(batch_size, )\n        \n        for j, row in enumerate(rows):\n            samples[j] = create_X(data[:, 0], last_index=row, n_steps=n_steps, step_length=step_length)\n            targets[j] = data[row - 1, 1]\n        yield samples, targets","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:38.239775Z","iopub.execute_input":"2022-04-06T03:45:38.240001Z","iopub.status.idle":"2022-04-06T03:45:38.247332Z","shell.execute_reply.started":"2022-04-06T03:45:38.239970Z","shell.execute_reply":"2022-04-06T03:45:38.245852Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Position of second (of 16) earthquake. Used to have a clean split\n# between train and validation\nbatch_size = 32\nsecond_earthquake = 50085877\nfloat_data[second_earthquake, 1]","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:38.248689Z","iopub.execute_input":"2022-04-06T03:45:38.248938Z","iopub.status.idle":"2022-04-06T03:45:38.257292Z","shell.execute_reply.started":"2022-04-06T03:45:38.248902Z","shell.execute_reply":"2022-04-06T03:45:38.256420Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Initialize generators\ntrain_gen = generator(float_data, batch_size=batch_size) # Use this for better score\n# train_gen = generator(float_data, batch_size=batch_size, min_index=second_earthquake + 1)\nvalid_gen = generator(float_data, batch_size=batch_size, max_index=second_earthquake)\n","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:38.258936Z","iopub.execute_input":"2022-04-06T03:45:38.259206Z","iopub.status.idle":"2022-04-06T03:45:38.264665Z","shell.execute_reply.started":"2022-04-06T03:45:38.259167Z","shell.execute_reply":"2022-04-06T03:45:38.263805Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import tensorflow as tf\n# Define model\nfrom keras.models import Sequential\nfrom keras.layers import Dense#, CuDNNGRU\nfrom tensorflow.compat.v1.keras.layers import CuDNNGRU\nfrom keras.callbacks import ModelCheckpoint","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:38.266184Z","iopub.execute_input":"2022-04-06T03:45:38.266466Z","iopub.status.idle":"2022-04-06T03:45:39.569411Z","shell.execute_reply.started":"2022-04-06T03:45:38.266431Z","shell.execute_reply":"2022-04-06T03:45:39.568646Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"cb = [ModelCheckpoint(\"model.hdf5\", save_best_only=True, period=3)]\n","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:39.570588Z","iopub.execute_input":"2022-04-06T03:45:39.570848Z","iopub.status.idle":"2022-04-06T03:45:39.580557Z","shell.execute_reply.started":"2022-04-06T03:45:39.570814Z","shell.execute_reply":"2022-04-06T03:45:39.577967Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model = Sequential()\nmodel.add(CuDNNGRU(48, input_shape=(None, n_features)))\nmodel.add(Dense(10, activation='relu'))\nmodel.add(Dense(1))\n\nmodel.summary()","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:45:39.582873Z","iopub.execute_input":"2022-04-06T03:45:39.583125Z","iopub.status.idle":"2022-04-06T03:45:42.155725Z","shell.execute_reply.started":"2022-04-06T03:45:39.583089Z","shell.execute_reply":"2022-04-06T03:45:42.155036Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from tensorflow.keras.optimizers import Adam\n# Compile and fit model\nmodel.compile(optimizer=Adam(lr=0.0005), loss=\"mae\")\n\nhistory = model.fit(train_gen,\n                              steps_per_epoch=1000,\n                              epochs=30,\n                              verbose=0,\n                              callbacks=cb,\n                              validation_data=valid_gen,\n                              validation_steps=200)","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:52:22.229880Z","iopub.execute_input":"2022-04-06T03:52:22.230343Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Visualize accuracies\nimport matplotlib.pyplot as plt\ndef perf_plot(history, what = 'loss'):\n    x = history.history[what]\n    val_x = history.history['val_' + what]\n    epochs = np.asarray(history.epoch) + 1\n    \n    plt.plot(epochs, x, 'bo', label = \"Training \" + what)\n    plt.plot(epochs, val_x, 'b', label = \"Validation \" + what)\n    plt.title(\"Training and validation \" + what)\n    plt.xlabel(\"Epochs\")\n    plt.legend()\n    plt.show()\n    return None","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:50:46.008044Z","iopub.execute_input":"2022-04-06T03:50:46.008302Z","iopub.status.idle":"2022-04-06T03:50:46.015202Z","shell.execute_reply.started":"2022-04-06T03:50:46.008272Z","shell.execute_reply":"2022-04-06T03:50:46.014203Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"perf_plot(history)","metadata":{"execution":{"iopub.status.busy":"2022-04-06T03:50:49.772942Z","iopub.execute_input":"2022-04-06T03:50:49.773881Z","iopub.status.idle":"2022-04-06T03:50:50.006715Z","shell.execute_reply.started":"2022-04-06T03:50:49.773829Z","shell.execute_reply":"2022-04-06T03:50:50.006035Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}