In this Notebook, I had implemented the anomaly detection model to predict whether the energy usage in a building is anomalous or not. Significance of this project lies in the fact that the anomalous usage implies energy wastage and we can prevent this by identifying such instances. This will help to save a lot of energy. In this notebook I had just used the train.csv data to build a model.
\n
TABLE OF CONTENTS
\n \n* [1. Importing Libraries](#1)\n \n* [2. Exploratory Data Analysis and Feature Extraction](#2)\n \n* [3. Equalizing Data and Preprocessing](#3)\n \n* [4. Model Development and Predictions](#4) \n \n* [5. The End](#5) \n\n
\n#
IMPORTING LIBRARIES
","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sb\nfrom tqdm import tqdm, trange\nfrom datetime import datetime\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.preprocessing import LabelEncoder, StandardScaler\nfrom sklearn import metrics\nfrom xgboost import XGBClassifier\nfrom sklearn.linear_model import LogisticRegression\n\nimport warnings\nwarnings.filterwarnings('ignore')","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:19.208112Z","iopub.execute_input":"2022-05-21T17:41:19.209003Z","iopub.status.idle":"2022-05-21T17:41:20.782172Z","shell.execute_reply.started":"2022-05-21T17:41:19.208863Z","shell.execute_reply":"2022-05-21T17:41:20.781327Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
\n#
Exploratory Data Analysis and Feature Extraction
\n
By using the timestamp we will develop some features and side by side make some visualization to predict the relation ship between different variables.
","metadata":{}},{"cell_type":"code","source":"train = pd.read_csv('../input/energy-anomaly-detection/train.csv')\ntrain.head()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:24.383280Z","iopub.execute_input":"2022-05-21T17:41:24.383542Z","iopub.status.idle":"2022-05-21T17:41:25.857382Z","shell.execute_reply.started":"2022-05-21T17:41:24.383514Z","shell.execute_reply":"2022-05-21T17:41:25.856405Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
Clearly there are null values in the meter_reading column.
","metadata":{}},{"cell_type":"code","source":"train.shape","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:31.178161Z","iopub.execute_input":"2022-05-21T17:41:31.178446Z","iopub.status.idle":"2022-05-21T17:41:31.185830Z","shell.execute_reply.started":"2022-05-21T17:41:31.178415Z","shell.execute_reply":"2022-05-21T17:41:31.184721Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train.isnull().sum()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:31.397960Z","iopub.execute_input":"2022-05-21T17:41:31.398645Z","iopub.status.idle":"2022-05-21T17:41:31.602558Z","shell.execute_reply.started":"2022-05-21T17:41:31.398590Z","shell.execute_reply":"2022-05-21T17:41:31.601873Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"k = train.isnull().sum().sum()\nplt.pie([k, train.shape[0]-k], labels=['Null', 'Filled']\\\n , colors=sb.color_palette('bright'), autopct='%.0f%%')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:31.677822Z","iopub.execute_input":"2022-05-21T17:41:31.678642Z","iopub.status.idle":"2022-05-21T17:41:32.001983Z","shell.execute_reply.started":"2022-05-21T17:41:31.678590Z","shell.execute_reply":"2022-05-21T17:41:32.000997Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\n
So, around 6% of the total entry's have null value in the meter_reading column.
","metadata":{}},{"cell_type":"code","source":"plt.subplots(figsize=(15,5))\nplt.subplot(1,2,1)\nsb.countplot(x=train['anomaly'])\nplt.subplot(1,2,2)\ndata = train['anomaly'].value_counts()\nplt.pie(data, labels=['Normal', 'Anomalous'], colors=sb.color_palette('bright'), autopct='%.0f%%')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:38.688729Z","iopub.execute_input":"2022-05-21T17:41:38.689032Z","iopub.status.idle":"2022-05-21T17:41:39.115673Z","shell.execute_reply.started":"2022-05-21T17:41:38.688998Z","shell.execute_reply":"2022-05-21T17:41:39.114978Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
As expected generally in a anomaly detection problem data imbalance is an obvious case because anomalies as rare cases which do not occur oftenly else there must be some other default.
","metadata":{}},{"cell_type":"code","source":"mean_reading = train.groupby('building_id').mean()['meter_reading']","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:44.738744Z","iopub.execute_input":"2022-05-21T17:41:44.739388Z","iopub.status.idle":"2022-05-21T17:41:44.832288Z","shell.execute_reply.started":"2022-05-21T17:41:44.739334Z","shell.execute_reply":"2022-05-21T17:41:44.831052Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(20,10))\nmean_reading.plot.bar()\nplt.xticks([])\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:44.968540Z","iopub.execute_input":"2022-05-21T17:41:44.969383Z","iopub.status.idle":"2022-05-21T17:41:45.845021Z","shell.execute_reply.started":"2022-05-21T17:41:44.969340Z","shell.execute_reply":"2022-05-21T17:41:45.844316Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"mean_reading.idxmax()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:45.846341Z","iopub.execute_input":"2022-05-21T17:41:45.846732Z","iopub.status.idle":"2022-05-21T17:41:45.852771Z","shell.execute_reply.started":"2022-05-21T17:41:45.846702Z","shell.execute_reply":"2022-05-21T17:41:45.851859Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
- The building with id equal to
801 has mean meter_reading too high as compare to other buildings seems like a skyscrapper to me😆 and some buildings seems to have very low mean meter_reading. - \n\nThis mean statistics is for throughout the year. We will see further in this notebook is there some variation in the mean usage during different days/month.
","metadata":{}},{"cell_type":"code","source":"def add_date_time(data):\n new = data['timestamp'].str.split(' ', n = 1, expand=True)\n data['date'] = new[0]\n data['time'] = new[1]\n data['time'] = data['time'].str[:2].astype(int)\n data = data.drop(['timestamp'], axis=1)\n\n return data\n\ntrain = add_date_time(train)\ntrain.head()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:41:56.089780Z","iopub.execute_input":"2022-05-21T17:41:56.090165Z","iopub.status.idle":"2022-05-21T17:42:01.945735Z","shell.execute_reply.started":"2022-05-21T17:41:56.090125Z","shell.execute_reply":"2022-05-21T17:42:01.944400Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(15,5))\ntrain.groupby('time').mean()['meter_reading'].plot.bar()\nplt.title('Mean usage hourly.', fontsize=20)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:01.947878Z","iopub.execute_input":"2022-05-21T17:42:01.948435Z","iopub.status.idle":"2022-05-21T17:42:02.460032Z","shell.execute_reply.started":"2022-05-21T17:42:01.948384Z","shell.execute_reply":"2022-05-21T17:42:02.458982Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def add_month_day(data):\n data['month'] = data['date'].str[5:7].astype(int)\n data['day'] = data['date'].str[-2:].astype(int)\n data = data.drop(['date'], axis=1)\n\n return data\n\ntrain = add_month_day(train)\ntrain.head()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:02.461443Z","iopub.execute_input":"2022-05-21T17:42:02.461768Z","iopub.status.idle":"2022-05-21T17:42:04.618843Z","shell.execute_reply.started":"2022-05-21T17:42:02.461727Z","shell.execute_reply":"2022-05-21T17:42:04.617817Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(15,5))\ntrain.groupby('month').mean()['meter_reading'].plot.bar()\nplt.title('Mean usage monthly.', fontsize=20)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:04.620854Z","iopub.execute_input":"2022-05-21T17:42:04.621110Z","iopub.status.idle":"2022-05-21T17:42:04.968660Z","shell.execute_reply.started":"2022-05-21T17:42:04.621079Z","shell.execute_reply":"2022-05-21T17:42:04.967598Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(15,5))\ntrain.groupby('day').mean()['meter_reading'].plot.bar()\nplt.title('Mean usage daily.', fontsize=20)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:04.969876Z","iopub.execute_input":"2022-05-21T17:42:04.970112Z","iopub.status.idle":"2022-05-21T17:42:05.419225Z","shell.execute_reply.started":"2022-05-21T17:42:04.970084Z","shell.execute_reply":"2022-05-21T17:42:05.418352Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train.head()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:05.421045Z","iopub.execute_input":"2022-05-21T17:42:05.421668Z","iopub.status.idle":"2022-05-21T17:42:05.434828Z","shell.execute_reply.started":"2022-05-21T17:42:05.421603Z","shell.execute_reply":"2022-05-21T17:42:05.433953Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Feature denoting whether it was a weekday or weekend.\ndef weekend_or_weekday(year,month,day):\n \n d = datetime(year,month,day)\n if d.weekday()>4:\n return 1\n else:\n return 0\n\n\ntrain['weekend'] = train.apply(lambda x:\\\n weekend_or_weekday(int(2016),int(x['month']),int(x['day'])),axis=1)\ntrain.head()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:05.436589Z","iopub.execute_input":"2022-05-21T17:42:05.436998Z","iopub.status.idle":"2022-05-21T17:42:40.641200Z","shell.execute_reply.started":"2022-05-21T17:42:05.436961Z","shell.execute_reply":"2022-05-21T17:42:40.640345Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train['weekend'].value_counts()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:40.642759Z","iopub.execute_input":"2022-05-21T17:42:40.643240Z","iopub.status.idle":"2022-05-21T17:42:40.662350Z","shell.execute_reply.started":"2022-05-21T17:42:40.643194Z","shell.execute_reply":"2022-05-21T17:42:40.661712Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data = train.groupby('weekend')['meter_reading'].mean()\nplt.pie(data, labels=['Weekday Mean Usage', 'Weekend Mean Usage']\\\n , colors=sb.color_palette('bright'))\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:40.663486Z","iopub.execute_input":"2022-05-21T17:42:40.663826Z","iopub.status.idle":"2022-05-21T17:42:40.791250Z","shell.execute_reply.started":"2022-05-21T17:42:40.663791Z","shell.execute_reply":"2022-05-21T17:42:40.790354Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
This means that the mean power usage of th buildings are higher during weekdays than the weekends which seems like obvious because offices, schools, colleges remains shut on the weekends.
","metadata":{}},{"cell_type":"code","source":"def impute_nulls(data):\n mean_reading = data.groupby('building_id').mean()['meter_reading']\n building_id = mean_reading.index\n values = mean_reading.values\n \n for i, idx in enumerate(building_id):\n data[data['building_id']==idx] = data[data['building_id']==idx].fillna(values[i]) \n \n return data\n\ntrain = impute_nulls(train)","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:40.794616Z","iopub.execute_input":"2022-05-21T17:42:40.795334Z","iopub.status.idle":"2022-05-21T17:42:56.315891Z","shell.execute_reply.started":"2022-05-21T17:42:40.795280Z","shell.execute_reply":"2022-05-21T17:42:56.315126Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
Although from the above charts we can observe that the mean usage of all the building is approximately same throughout the year (monthly, daily) and is around 175 kwh approx. but we cannot impute the missing values with this as some buildings have very high usage stats and some buildings have readings even close to zero.
","metadata":{}},{"cell_type":"markdown","source":"
\n#
Equalizing Data and Preprocessing
\n
As there is a high imbalance in the data so, if we won't equalize the positive and the negative examples in the data then the model will have a hard time generalizing the positive examples which is the main goal here.
","metadata":{}},{"cell_type":"code","source":"neg = train[train['anomaly'] == 0]\npos = train[train['anomaly'] == 1]\n\nprint(neg.shape, pos.shape)\nnegs1 = neg.sample(n = 37296, random_state=10)\nnegs2 = neg.sample(n = 37296, random_state=20)\ndf_eq = pd.concat([negs1, pos, negs2, pos], axis=0)\nprint(df_eq.shape)","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:56.317014Z","iopub.execute_input":"2022-05-21T17:42:56.317368Z","iopub.status.idle":"2022-05-21T17:42:56.570097Z","shell.execute_reply.started":"2022-05-21T17:42:56.317337Z","shell.execute_reply":"2022-05-21T17:42:56.568828Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
Yes, by using this method our model won't be able to see all the data even once but this will helps the model to generalize well on the positive as well as negative examples.
","metadata":{}},{"cell_type":"code","source":"features = df_eq.drop(['building_id', 'anomaly'], axis=1)\ntarget = df_eq['anomaly']","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:56.571461Z","iopub.execute_input":"2022-05-21T17:42:56.572274Z","iopub.status.idle":"2022-05-21T17:42:56.585188Z","shell.execute_reply.started":"2022-05-21T17:42:56.572231Z","shell.execute_reply":"2022-05-21T17:42:56.584468Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train, X_val, Y_train, Y_val = train_test_split(features, target, test_size = 0.01, random_state=10)\nprint(X_train.shape, X_val.shape)\n\nscaler = StandardScaler()\nX_train = scaler.fit_transform(X_train)\nX_val = scaler.transform(X_val)","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:56.586308Z","iopub.execute_input":"2022-05-21T17:42:56.586939Z","iopub.status.idle":"2022-05-21T17:42:56.633525Z","shell.execute_reply.started":"2022-05-21T17:42:56.586903Z","shell.execute_reply":"2022-05-21T17:42:56.632405Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
\n#
Model Development and Predictions
\n
Training data is ready now so, let us use some state of the art model to get a baseline on which we can build up by using the features data. Same set of processing needs to be done in the same order as we did in the training data so, that we don't do any mistake.
","metadata":{}},{"cell_type":"code","source":"models = [LogisticRegression(), XGBClassifier()]\n\nfor i in range(2):\n models[i].fit(X_train, Y_train)\n\n print(f'{models[i]} : ')\n print('Training Accuracy : ', metrics.roc_auc_score(Y_train, models[i].predict_proba(X_train)[:,1]))\n print('Validation Accuracy : ', metrics.roc_auc_score(Y_val, models[i].predict_proba(X_val)[:,1]))\n print()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:42:56.635015Z","iopub.execute_input":"2022-05-21T17:42:56.635280Z","iopub.status.idle":"2022-05-21T17:43:02.314316Z","shell.execute_reply.started":"2022-05-21T17:42:56.635248Z","shell.execute_reply":"2022-05-21T17:43:02.313460Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
- Logistic regressions accuracy is too low but
XGBClassifier is performing really well without any overfitting it give 94% validation accuracy. - \nAt this point you may have one doubt in your mind that why I am predicting probabilities instead of binary predictions then please pay a visit to this discussion do let me know your thoughts about it.
","metadata":{}},{"cell_type":"code","source":"test = pd.read_csv('../input/energy-anomaly-detection/test.csv', index_col=0)\ntest = add_date_time(test)\ntest = add_month_day(test)\ntest['weekend'] = test.apply(lambda x:\\\n weekend_or_weekday(int(2016),int(x['month']),int(x['day'])),axis=1)\ntest.head()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:43:02.316584Z","iopub.execute_input":"2022-05-21T17:43:02.316851Z","iopub.status.idle":"2022-05-21T17:43:51.473655Z","shell.execute_reply.started":"2022-05-21T17:43:02.316821Z","shell.execute_reply":"2022-05-21T17:43:51.472647Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test.shape","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:43:51.475084Z","iopub.execute_input":"2022-05-21T17:43:51.475396Z","iopub.status.idle":"2022-05-21T17:43:51.482549Z","shell.execute_reply.started":"2022-05-21T17:43:51.475353Z","shell.execute_reply":"2022-05-21T17:43:51.481559Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test.isnull().sum()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:43:51.483969Z","iopub.execute_input":"2022-05-21T17:43:51.484263Z","iopub.status.idle":"2022-05-21T17:43:51.516274Z","shell.execute_reply.started":"2022-05-21T17:43:51.484224Z","shell.execute_reply":"2022-05-21T17:43:51.515411Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"k = test.isnull().sum().sum()\nplt.pie([k, test.shape[0]-k], labels=['Null', 'Filled']\\\n , colors=sb.color_palette('bright'), autopct='%.0f%%')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:43:51.517740Z","iopub.execute_input":"2022-05-21T17:43:51.518010Z","iopub.status.idle":"2022-05-21T17:43:51.748686Z","shell.execute_reply.started":"2022-05-21T17:43:51.517977Z","shell.execute_reply":"2022-05-21T17:43:51.747718Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
So, 5% of the entries in the meter_reading of the test data is null which is similar to the train data which has approximately 6% of the meter_reading as null.
","metadata":{}},{"cell_type":"code","source":"test = impute_nulls(test)\ntest.isnull().sum().sum()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:43:51.750295Z","iopub.execute_input":"2022-05-21T17:43:51.751552Z","iopub.status.idle":"2022-05-21T17:44:06.298410Z","shell.execute_reply.started":"2022-05-21T17:43:51.751487Z","shell.execute_reply":"2022-05-21T17:44:06.297560Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test = test.drop(['building_id'], axis=1)\ntest = scaler.transform(test)","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:44:06.299793Z","iopub.execute_input":"2022-05-21T17:44:06.300118Z","iopub.status.idle":"2022-05-21T17:44:06.409416Z","shell.execute_reply.started":"2022-05-21T17:44:06.300073Z","shell.execute_reply":"2022-05-21T17:44:06.408317Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"predictions = models[1].predict_proba(test)[:,1]","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:44:06.410993Z","iopub.execute_input":"2022-05-21T17:44:06.411340Z","iopub.status.idle":"2022-05-21T17:44:08.157468Z","shell.execute_reply.started":"2022-05-21T17:44:06.411298Z","shell.execute_reply":"2022-05-21T17:44:08.156799Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ss = pd.read_csv('../input/energy-anomaly-detection/sample_submission.csv')\nss['anomaly'] = predictions\nss.to_csv('Submission.csv', index=False)\nss.head()","metadata":{"execution":{"iopub.status.busy":"2022-05-21T17:44:08.158762Z","iopub.execute_input":"2022-05-21T17:44:08.159192Z","iopub.status.idle":"2022-05-21T17:44:13.607115Z","shell.execute_reply.started":"2022-05-21T17:44:08.159143Z","shell.execute_reply":"2022-05-21T17:44:13.606194Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"
\n#
The End
\n
- Score evaluaton using the binary prediction and probability prediction gives different scores on the same set of predictions do visit the above discussion to know why this happens.
- I am working on another notebook as well in which I will develop a model by using the
train_features.csv data to get better results.
","metadata":{}}]}