{"cells":[{"metadata":{},"cell_type":"markdown","source":"Hey guys!\n\nAfter the competition closed, I want to share my idea of filling missing with you.\n(The original train data have too mcuh missing...)\n\nThe main idea is \"time-series similarity between meters in single site\".\n(for example, meter_A's pattern is very similar to meter_B in site1, and I do believe that some of you have noticed that)\nBy well utilizing such similarity, we can find out those missing data!"},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install kneed","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":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 lightgbm as lgb\nfrom pathlib import Path\n\nimport matplotlib.pyplot as plt\nimport plotly.graph_objs as go\nimport plotly\nfrom plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot\nimport cufflinks as cf\ncf.set_config_file(offline=True)\n\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\nfrom statsmodels.tsa.seasonal import seasonal_decompose\n\nfrom kneed import KneeLocator\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 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# Any results you write to the current directory are saved as output.","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def MinMaxScaler(data):\n    return (data-np.min(data))/(np.max(data)-np.min(data))\n\ndef Kmeans_clustering(df, clusterNum, max_iter, n_jobs):\n    scaler = StandardScaler()\n    scaler.fit(df)\n    df_std = pd.DataFrame(data=scaler.transform(df), columns=df.columns, index=df.index)\n    km_model = KMeans(n_clusters=clusterNum, max_iter=max_iter, n_jobs=n_jobs, random_state=666)\n    km_model = km_model.fit(df_std)\n    clusterdf= pd.DataFrame(data=km_model.labels_, columns=['ClusterNo'])\n    clusterdf.index = df.index\n    return clusterdf\n\ndef Kmeans_bestClusterNum(df, range_min, range_max, max_iter, n_jobs):\n    silhouette_avgs = []\n    sum_of_squared_distances = []\n    \n    ks = range(range_min,range_max+1)\n    for k in ks:\n        kmeans_fit = KMeans(n_clusters = k, n_jobs=n_jobs, max_iter=max_iter, random_state=666).fit(df)\n        cluster_labels = kmeans_fit.labels_\n        sum_of_squared_distances.append(kmeans_fit.inertia_)\n        \n    kn = KneeLocator(list(ks), sum_of_squared_distances, S=1.0, curve='convex', direction='decreasing')  \n    plt.xlabel('k')\n    plt.ylabel('sum_of_squared_distances')\n    plt.title('The Elbow Method showing the optimal k')\n    plt.plot(ks, sum_of_squared_distances, 'bx-')\n    plt.vlines(kn.knee, plt.ylim()[0], plt.ylim()[1], linestyles='dashed')\n    print('Optimal clustering number:'+str(kn.knee))\n    print('----------------------------')    \n    \n    return kn.knee","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"%%time\nroot = Path('../input/ashrae-feather-format-for-fast-loading')\n\ntrain_df = pd.read_feather(root/'train.feather')\ntest_df = pd.read_feather(root/'test.feather')\nweather_train_df = pd.read_feather(root/'weather_train.feather')\nweather_test_df = pd.read_feather(root/'weather_test.feather')\nbuilding_meta_df = pd.read_feather(root/'building_metadata.feather')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#Let's just take site1 as example\ntrain_df = train_df.merge(building_meta_df, on='building_id')\ntrain_df = train_df[train_df['site_id']==1]\n\ntrain_df['merged_id'] = 'site'+train_df['site_id'].astype(str)+'_'\\\n                        +'bldg'+train_df['building_id'].astype(str)+'_'\\\n                        +'meter'+train_df['meter'].astype(str)\ntrain_df = train_df[['timestamp', 'site_id', 'building_id', 'meter' ,'meter_reading', 'merged_id']]\ntrain_df = train_df.set_index('timestamp')\n\ntrain_df = train_df.sort_values(['merged_id','timestamp'])\ntrain_df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_df_pivot = train_df.pivot_table(values='meter_reading', index='timestamp', columns='merged_id')\ntrain_df_pivot","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df_PM_temp = train_df_pivot.copy()\ndf_PM_temp = (df_PM_temp-df_PM_temp.mean())/df_PM_temp.std()\nSTL_decomp = seasonal_decompose(df_PM_temp.fillna(0), model='additive',freq=24*7,extrapolate_trend=1)\n\ndf_seasonal_temp=STL_decomp.seasonal.iloc[24*10:24*17,:]\ndf_seasonal_temp = df_seasonal_temp.T\n\ndf_seasonal_temp['ClusterNo'] = Kmeans_clustering(df=df_seasonal_temp, clusterNum=5, max_iter=100000, n_jobs=-1)\n\nfor ClusterNo in df_seasonal_temp['ClusterNo'].unique():\n    df_plot = df_seasonal_temp[df_seasonal_temp['ClusterNo']==ClusterNo].T.drop('ClusterNo')\n    print('ClusterNo: ' + str(ClusterNo))    \n    print('Amount of meters: ' + str(len(df_plot.T)))\n    df_plot.plot(figsize=(15,5),color='black',alpha=0.1,legend=False)\n    plt.show()\n    print('---------------------------------------------------------------------------------------------------')    ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df_PM_temp = train_df_pivot.copy()\ndf_PM_temp = (df_PM_temp-df_PM_temp.mean())/df_PM_temp.std()\nSTL_decomp = seasonal_decompose(df_PM_temp.fillna(0), model='additive',freq=24*7,extrapolate_trend=1)\n\ndf_trend_temp = STL_decomp.trend\n#df_trend_temp = df_trend_temp.loc['2016']\ndf_trend_temp = df_trend_temp.T\n\ndf_trend_temp['ClusterNo'] = Kmeans_clustering(df=df_trend_temp, clusterNum=5, max_iter=100000, n_jobs=-1)\n\nfor ClusterNo in df_trend_temp['ClusterNo'].unique():\n    df_plot = df_trend_temp[df_trend_temp['ClusterNo']==ClusterNo].T.drop('ClusterNo')\n    print('ClusterNo: ' + str(ClusterNo))    \n    print('Amount of meters: ' + str(len(df_plot.T)))\n    df_plot.plot(figsize=(15,5),color='black',alpha=0.1,legend=False,ylim=(-3,3))\n    plt.show()\n    print('---------------------------------------------------------------------------------------------------')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Take a look of \"site1_bldg105_meter0\", we can observe abnormal constant in some days."},{"metadata":{"trusted":true},"cell_type":"code","source":"#Graph of original data: \"site1_bldg105_meter0\"\ntrain_df_pivot['site1_bldg105_meter0'].iplot(kind='scatter')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#Let's clean these abnormal days by checking std_daily==0\n\nlist_powerMeter = list(train_df_pivot.columns)\n\nfor name_meter in list_powerMeter:\n    df_daily = train_df_pivot[name_meter].copy()\n    df_daily = df_daily.resample('D').std()\n\n    list_abnormalDate = list(df_daily[df_daily==0].index.strftime('%Y-%m-%d'))\n    train_df_pivot['Date'] = pd.to_datetime(train_df_pivot.index.date)\n\n    train_df_pivot.loc[train_df_pivot['Date'].isin(list_abnormalDate), name_meter] = np.nan\n    train_df_pivot = train_df_pivot.drop('Date', axis=1)  ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Successfully clean these abnormal days!"},{"metadata":{"trusted":true},"cell_type":"code","source":"#Graph of cleand data: \"site1_bldg105_meter0\"\ntrain_df_pivot['site1_bldg105_meter0'].iplot(kind='scatter')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Now, let's take \"site1_bldg105_meter0\" as example for fillling missing"},{"metadata":{"trusted":true},"cell_type":"code","source":"#Prepare data: \"site1_bldg105_meter0\"\nexample_df = train_df_pivot.copy()\nexample_df = example_df.dropna(axis=1) #Only take non-missing meter readings as inputs\n\nexample_df['elec_meas'] = train_df_pivot['site1_bldg105_meter0'].copy() #Take \"site1_bldg105_meter0\" as target\n\ntraindata = example_df[~example_df['elec_meas'].isna()]\ntestdata = example_df[example_df['elec_meas'].isna()]\n\ntrain_labels = traindata['elec_meas']\n\ntrain_features = traindata.drop('elec_meas', axis=1)\ntest_features = testdata.drop('elec_meas', axis=1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"LGB_model = lgb.LGBMRegressor()\nLGB_model.fit(train_features, train_labels)\n\ntestdata['elec_pred'] = LGB_model.predict(test_features)\n\n# Use the forest's predict method on the train data\nexample_df['elec_pred'] = LGB_model.predict(example_df.drop('elec_meas', axis=1))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Boom! Here's the meter reading after filling missing😎"},{"metadata":{"trusted":true},"cell_type":"code","source":"example_df[['elec_meas', 'elec_pred']].iplot(kind='scatter')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Hope you guys enjoy!"}],"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":4,"nbformat_minor":1}