{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"markdown","source":"![](http://retrofittingca.com/wp-content/uploads/2017/02/PROOFED-RetrofittingCA-WhyEarthquakesAreSoHardtoPredictJanuary312017-PIC.jpg)\n"},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","collapsed":true,"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":false},"cell_type":"markdown","source":"# How to handle memory"},{"metadata":{"_uuid":"598742f732017df0794ccada8e1c38eb6a3f924f"},"cell_type":"markdown","source":" <a id=\"top\"></a> <br>\n## Notebook  Content\n1. [Memory issues](#1)\n1. [EDA-briefly](#2)\n1. [Feature Engineering (to be continued) ](#3)\n"},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"168e1de59a0e54e2199ae95b86420201c7461ec0"},"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport os\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\nimport warnings\nwarnings.filterwarnings(\"ignore\")\nimport gc","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"2169376d8ccea60eb9f52f5baef42b95aafe62f6"},"cell_type":"markdown","source":"<a id=\"1\"></a> <br>\n#  1-Memory issues"},{"metadata":{"trusted":true,"_uuid":"116d12bfdb248def2ec52d21f2c110ecc5d5f1a4"},"cell_type":"code","source":"%%time\ntrain = pd.read_csv('../input/LANL-Earthquake-Prediction/train.csv', dtype={'acoustic_data': np.int16, 'time_to_failure': np.float64})","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d3cdd30469dfa55283e6e684dc5ce0161f279970"},"cell_type":"markdown","source":"First potential problem will be definately memory usage. WE have over 629 **Milion** rows of date, and as we can see even with pre-setting data-types it takes a couple of minutes. Right of the bat 10GB is taken. We can fight it with following approaches:\n\n*  **DASK** ---(\"A Dask DataFrame is a large parallel dataframe composed of many smaller Pandas dataframes, split along the index. These pandas dataframes may live on disk for larger-than-memory computing on a single machine, or on many different machines in a cluster. One Dask dataframe operation triggers many operations on the constituent Pandas dataframes.\") [Dask-doc]((http://docs.dask.org/en/latest/dataframe.html) **NUTSHELL:**  python library for parallel computing that can work on a single notebook or large cluster and will compute the operations (careful not all of them) very fast\n* ** Garbage Collector** (gc()) and deleting unused variables\n*  **Pre-setting data types** when importing (as above)\n* **Custom Function** to reduce memory (to be found in 2. EDA section)\n\n"},{"metadata":{"trusted":true,"_uuid":"793d04c4c05e58109d1e556de30a261ef7bc3fce"},"cell_type":"code","source":"import dask\nimport dask.dataframe as dd","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"df4ebcc3aeee487992dfd4dfb461625396af5b71"},"cell_type":"code","source":"%%time\ntrain1 = dd.read_csv(\"../input/LANL-Earthquake-Prediction/train.csv\",  dtype={'acoustic_data': np.int16, 'time_to_failure': np.float64})\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"c86b51517aef8943e0a6190632fbd917e6c2e696"},"cell_type":"markdown","source":"Lighting fast!. WE do have to use .compute() method every time, and not all pandas operations are covered :( "},{"metadata":{"trusted":true,"_uuid":"7b20127ebdcee0996baf5af457b56d21b9c110c5"},"cell_type":"code","source":"train.head()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"339e4e3585325c434ec0380d95b1ed241adc895a"},"cell_type":"markdown","source":"As already elaborated, given the accoustin data we ought to predict time  (in seconds) until the next earthquake takes place.** One peculiar thing** that will be important for feature engineering is that train set has a continous measurment, while test set has a **2624** Sets. Each with 150 000 rows (150 000 acoustic data) and we were supposed to predict the time until failure. It is intersting because these are only snippets of sime series of data and are (possibly) not even connected to each other. One implication of such framework is that we have to engineer our new variables on different (random) snippets of the accoustic_data series."},{"metadata":{"_uuid":"0a0d28865240b37e71ef3bc15585a505d8173cd4"},"cell_type":"markdown","source":"<a id=\"2\"></a> <br>\n#  2-EDA\nAlready some great EDA analysis was made, I will just show the most interesting (in my opinion) relationships that can be noticed."},{"metadata":{"trusted":true,"_uuid":"9c7ed1d722cc6a265b8c4994b52f7fc872c5de3f"},"cell_type":"code","source":"plot1=train1[0:150000].time_to_failure\nplot2=train1[0:150000].acoustic_data\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"fe6a55d1ee2de8a28086aaea0006ab89f42967f7"},"cell_type":"markdown","source":"Why 150 000? That is the size of one test example that we ought to predict. ANd why negative bvalues of time to failure? In order to \"approach\" the zero..."},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"f588e7a762d473f1925ee0da15a63790f5acb767"},"cell_type":"code","source":"%%time\nplt.figure(figsize=(20,5))\nplt.plot(-plot1.compute(), plot2.compute());\nplt.xlabel(\"- Quaketime\")\nplt.ylabel(\"Signal\")\nplt.title(\"PLOT 0\");","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"96aa2d221cd32f7e1b7148b7b1ea33a12894045e"},"cell_type":"code","source":"%%time\nplt.figure(figsize=(20,5))\nplt.plot(-train.time_to_failure.values[0:150000], train.acoustic_data.values[0:150000]);\nplt.xlabel(\"- Quaketime\")\nplt.ylabel(\"Signal\")\nplt.title(\"PLOT 1\");","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"4823b152c936789723e744ae67b3fbf9517931f4"},"cell_type":"markdown","source":"# Conclusion of first plots:\nBoth from PLOT 1 and 0 we can see that there is huge change in signal right before we start going towards towards null. One should notice different values on x and y axis of PLOT 1 and 0. Plot 1 is more accurate in a sense that actually these values are to be found in the original dataframe. Dask, seems to depict values differently, but still observation is the same. Change in Signal value influences Quaketime"},{"metadata":{"_uuid":"3c5d897391d0702fd6b9321407a425e226d1272a"},"cell_type":"markdown","source":"# **Dask or Pandas**\nFrom the simple analysis above we can already see trade-offs in these two different approaches. It is different, it seems that values are differently stored, it is a different object type after all and a more efficient way to acces these values has to be out there. And run-time is horrible. As I said somethings are not implemented and user has to think about a work around. One can think about Dask as new pandas, with a bit more complicated synthax but upside regarding memory.\n\n"},{"metadata":{"_uuid":"acaac9b7344d7ec1985653ef6d6653cb14b0b66f","trusted":true},"cell_type":"code","source":"from IPython.display import Image\nImage(\"../input/example/example.png\")","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"1af07175cd58ae436daa791d26dbf080941723e4"},"cell_type":"markdown","source":"Ok so the learning curve is a bit steep but not **REALLY**.\n[Dask Tutorial](https://github.com/dask/dask-tutorial/blob/master/01_dask.delayed.ipynb) can help us get to date, and I suspect that winning kernels will most likely use dask."},{"metadata":{"_uuid":"bfc604c186051ed32926e84b002aa8e3a05d04e3"},"cell_type":"markdown","source":"# Conclusion regarding **Dask**\nDepending on the appetite and analysis that will follow in next months one should/should not use Dask. I will explicitly forgo it on this public kernel in order for other users to fork it maybe attach one owns analysis without thinking how should one modify it as a Dask data object. In that spirit one should definately do:"},{"metadata":{"_uuid":"0f1302c162cca321e79c6028c0e09779cb62738a","trusted":true},"cell_type":"code","source":"del train1\ngc.collect()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d1ead2d19df98a4cd7db79d2beccbd8cfb0bea64"},"cell_type":"markdown","source":"**Custom Function** Having only 2 columns, it isnt a deal breaker, but it helps. Functionality: iterate through all the columns of a dataframe and modify the datatype to reduce memory usage."},{"metadata":{"trusted":true,"_uuid":"6b3caef6382bd0008a8761bb8b32a6adeff5ef54"},"cell_type":"code","source":"def reduce_mem_usage(df):\n\n    for col in df.columns:\n        col_type = df[col].dtype.name\n\n        if col_type not in ['object', 'category', 'datetime64[ns, UTC]']:\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    return df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"66256f6843ea148a67d2d22decb668922c754067"},"cell_type":"code","source":"train=reduce_mem_usage(train)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"3452f85ef2ee8b4c9a28a0cba31b98b169b13836"},"cell_type":"code","source":"gc.collect()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"563640c95942b820b8a56331068d4b0f3fe904b3"},"cell_type":"markdown","source":"Let us proceede with further interesting **EDA**"},{"metadata":{"_kg_hide-input":false,"trusted":true,"_uuid":"adb103c148e016dc4b80c8524ef3b30b593b5041"},"cell_type":"code","source":"fig, ax1 = plt.subplots(figsize=(16, 8))\nplt.title(\"Trends of acoustic_data and time_to_failure. 150 000 rows\")\nplt.plot(train['acoustic_data'].values[:15000000], color='y')\nax1.set_ylabel('acoustic_data', color='y')\nplt.legend(['acoustic_data'])\nax2 = ax1.twinx()\nplt.plot(train['time_to_failure'].values[:15000000], color='r')\nax2.set_ylabel('time_to_failure', color='r')\nplt.legend(['time_to_failure'], loc=(0.875, 0.9))\nplt.grid(False)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"86b951c95c5707f5fdfdf03a228e38a975c481ba"},"cell_type":"markdown","source":"This time around we plotted the independently but on the same graph and the conclusion is the same. Right before the problem happens signal increases. Ok what wcan we now say regarding the univariate distribution of the columns themselves?"},{"metadata":{"trusted":true,"_uuid":"5b22ff113fc9a9f4629b4a0df86115a8b3dc8fa5"},"cell_type":"code","source":"fig, ax = plt.subplots(1,2, figsize=(20,5))\nsns.distplot(train['acoustic_data'].values[:15000000], ax=ax[0], color=\"Yellow\", bins=100)\nax[0].set_xlabel(\"acoustic_data\")\nax[0].set_ylabel(\"Density\")\nax[0].set_title(\"acoustic_data distribution\")\n\n\nsns.distplot(train['time_to_failure'].values[:15000000], ax=ax[1], color=\"Red\", bins=100)\nax[1].set_xlabel(\"time_to_failure\")\nax[1].set_ylabel(\"Density\")\nax[1].set_title(\"time_to_failure distribution\");","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"cdf35e7ba3d6f53f6069d380c7378ca1ccfd39f2"},"cell_type":"markdown","source":"Very intersting (again for the first 150 000 observations) it seems that (yellow=accoustic_data) follows normal distribution with some outliers. On the other hand time to failure takes on 2(?) very distinct values..."},{"metadata":{"_uuid":"ab17bf2805dd5738d5f98b5b0f87d97ca5f3290c"},"cell_type":"markdown","source":"# **Implications of EDA**\n\nThis EDA helps us actually with **feature engineering**. If we look at the **first** graph(yellow) and the same thing we could have noticed before is that acoustic_data has a huge spike in value just before the earthquake is about the happen. So a new variable could be a column that assigns a big value (say 100) when this huge values occurs (say 3 points of standard deviation) and 0 otherwise. We hope that te algorithm will pick up this indicator and it will help with predictions. We can test this later on with Variable Importance using lgbm.\n\n**Second** graph gives us another theory. We see that we have multiple (atleast 2) distinct centers where time_to_failure values are centered. That  also wasnt hard to spot from the graph before that. We see that the red line decreases linearly towards zero, than it has a sharp jump towards another count-down time until the next earthquake. Hence we can see that there are atleast two jumps/differences between time_to_failure values. First one is the slow linear differences as we approach the 0, and the next one is the steep jump. Implication is easy, create a feature that represents these categorical differences in time_to_failure/target variable."},{"metadata":{"_uuid":"538c9d79512a0b84473183e0c14f187aa5af5690"},"cell_type":"markdown","source":"<a id=\"3\"></a> <br>\n#  3-Feature Engineering"},{"metadata":{"_uuid":"3ccc3c9ac0f057980ac99fb600a302e218222937"},"cell_type":"markdown","source":"So with help of the initial analysis we know approximately what we sould do. Let us first tackle the independent variable \"acoustic_data\" that is what are all of the features that can be engineered from this variable."},{"metadata":{"_uuid":"ecbca0ae00ddeaa8362a8a55907a2a538ed85f88"},"cell_type":"markdown","source":"**a) Sharp-rise indicator**\nAs already discussed we want an indicator that tells explicitly about the sharp rise in value"},{"metadata":{"trusted":true,"_uuid":"e8573a5da7eb1517c95f01790e1fadd86df70031"},"cell_type":"code","source":"min_ = train.acoustic_data.mean() - 3 * train.acoustic_data.std()\nmax_ = train.acoustic_data.mean() + 3 * train.acoustic_data.std() ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"8902700b5d01d24f7f403720db3a75a25bc6c1fd"},"cell_type":"code","source":"#train[\"sharp_rise\"] = np.where((train.acoustic_data >= min_) & (train.acoustic_data <= max_), 100, 0)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"0cdd0368f28c4d11c414a64f67de394fee6f7521"},"cell_type":"markdown","source":"# **PART 1**\nIn order to showcase some things I did burn a lot of memory (I had to read whole data set without specifying nrows= etc...) **Part 2** of the same kernel should focus more on the actual results, rather exploration and memory discussion. It can be found right away at [Part 2](https://www.kaggle.com/zikazika/initial-approach-earthquake-prediction-part-2) I am sorry for the inconvenience."},{"metadata":{"_uuid":"a698834f8eab2401eaf6380245c1bd30db043f32"},"cell_type":"markdown","source":"**Upvotes are appreciated if you benefited :)**"},{"metadata":{"trusted":true,"_uuid":"b8ae595417a5c0faa4cdff27e2783d377db6ec19"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"afd7a125570890effe57db59ab58b9bae9dd323a"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"3389033543a6dc64a57486e272f4c387edfe7443"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"75b0c961e2fb2fb875d50eaae1f18bd4d4c95265"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"86bc666aeafbda0cd4c1cf78d5ba52afb0d5d3c8"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"431446714527b28d95b0341dc4fb2abced705d12"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e2e4b2d62c11370f4eb67c597c1cd7037d5bdead"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"2cd2f660cb1b924a276a78de3f98cc157a01492b"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d7138e8b0732c09bc9f2a68838528fd86c55f611"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"9b80d009e87e1d85e2c960aa789d21845a62bc30"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"be405a3247c658a0fe582d2379af78f4813d1eb7"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"6cf1f0f83648401815b7b1ead5b52f916c87e4ef"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4b38d412be984e243d26674da4ab0aa98ca071c0"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"7eab5a69d00ee87a0b080cdc0166aa2394fb2b66"},"cell_type":"markdown","source":""}],"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}