{"cells":[{"metadata":{"trusted":true},"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport seaborn as sns\nfrom sklearn.preprocessing import LabelEncoder\nimport datetime as dt \n%matplotlib inline\nimport matplotlib.pyplot as plt\nfrom tqdm import tqdm\nimport lightgbm as lgb\nimport datetime as dt\n\npd.options.display.float_format = '{:.4f}'.format","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"root = \"../input/ashrae-energy-prediction/\"\n\n#read data from csv files\ntrain = pd.read_csv(root + 'train.csv')\ntest = pd.read_csv(root + 'test.csv')\nw_train = pd.read_csv(root + 'weather_train.csv', index_col=False)\nw_test = pd.read_csv(root + 'weather_test.csv', index_col=False)\nbuilding = pd.read_csv(root + 'building_metadata.csv')\n\n#change the categorical primary use column to numerical with label encoder\nle = LabelEncoder()\nbuilding.primary_use = le.fit_transform(building.primary_use)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## Memory Reducer\n# :df pandas dataframe to reduce size             # type: pd.DataFrame()\n# :verbose                                        # type: bool\ndef reduce_mem_usage(df, verbose=True):\n    numerics = ['int16', 'int32', 'int64', 'float16', 'float32', 'float64']\n    start_mem = df.memory_usage().sum() / 1024**2    \n    for col in df.columns:\n        col_type = df[col].dtypes\n        if col_type in numerics:\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    end_mem = df.memory_usage().sum() / 1024**2\n    if verbose: print('Mem. usage decreased to {:5.2f} Mb ({:.1f}% reduction)'.format(end_mem, 100 * (start_mem - end_mem) / start_mem))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"reduce_mem_usage(train, verbose=True)\nreduce_mem_usage(test, verbose=True)\nreduce_mem_usage(w_train, verbose=True)\nreduce_mem_usage(building, verbose=True)\nreduce_mem_usage(w_test, verbose=True)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Merge"},{"metadata":{"trusted":true},"cell_type":"code","source":"train = train.merge(building, on='building_id', how='left')\ntrain = train.merge(w_train, on=['site_id', 'timestamp'], how='left')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test = test.merge(building, on='building_id', how='left')\ntest = test.merge(w_test, on=['site_id', 'timestamp'], how='left')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# We will need this row for our submission, \n#but it is not in the training dataset.\n#We want to drop it for our analysisso we will assign it a variable to recall later \nrow_id = test.row_id","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test.drop(columns=['row_id'], axis=1, inplace=True)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# EDA"},{"metadata":{"trusted":true},"cell_type":"code","source":"train.describe()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test.describe()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The dataset contains measurement error. We will need to look for ouliers produced by error. Over 9% of the meter readings are zero, but a building in use should not have a meter reading of zero. Either the reading is an error, or the meter was not read at that particular time. In the latter case, we will need to gather the weather data between non-zero readings to arrive at an accurate assessment of the readings."},{"metadata":{"trusted":true},"cell_type":"code","source":"print(\"Percent equal to zero:\", (train[train.meter_reading == 0].shape[0] / len(train.meter_reading))*100)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The meter_reading column records the difference energy use, not the aggregate as domestic electric meter does, where any reading must be taken as the difference between the prior reading to arrive at interval usage. We will take a look at buiding_id 2 to get a sense of the readings."},{"metadata":{"trusted":true},"cell_type":"code","source":"# An example of meter readings for building_id 2\ntrain[(train.building_id==2) & train.meter_reading>0].iloc[:100]","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"While considering the outliers, we will take a look at the percentage of meter readings above the mean."},{"metadata":{"trusted":true},"cell_type":"code","source":"print(\"Percent above mean:\", (train[train.meter_reading > train.meter_reading.mean()].shape[0] / train.shape[0]) * 100)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train.groupby('building_id').meter_reading.mean().sort_values()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"A better assessment of outliers would be to take the meter reading as a proportion of the building size, and we will tackle this in the Beta version."},{"metadata":{},"cell_type":"markdown","source":"The following shows the percent of nan's in each column for the train and test dataframes. We will drop those with a percentage."},{"metadata":{"trusted":true},"cell_type":"code","source":"train.isnull().sum() * 100 / len(test)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test.isnull().sum() * 100 / len(test)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Visualize data\n##### just a brief word on a few variables\nThe scatter matrices that follow are samples of 5000 data points. The target variable is the meter reading column. This version includes all meter types (elec, chilledwater). In the Beta version, we will explore the correlations within each meter type."},{"metadata":{},"cell_type":"markdown","source":"Dew temperature is the air temperature at which water in the air condenses. The higher the dew temperature the more \"muggy\" it feels. The data approaches a normal distribution (see, https://www.csemag.com/articles/controlling-dew-point/)"},{"metadata":{},"cell_type":"markdown","source":"Sea-level pressure lowers as altitude rises. Because the air is thinner at higher altitudes, it is less effective at heating at building. It takes more energy to heat and cool a building at a higher altitude than it does at sea level. For sea-level pressure, see https://www.britannica.com/science/atmospheric-pressure; https://w1.weather.gov/glossary/index.php?word=sea+level+pressure. For altitude and energy costs, see https://www.achrnews.com/articles/110408-selecting-heating-cooling-units-for-high-altitude-homes; http://www.comairrotron.com/solving-high-altitude-cooling-problems"},{"metadata":{"trusted":true},"cell_type":"code","source":"#create a sample dataframe that filters our values greater than 5000 and \n# not equal to 0, so that our vizualizations are to scale\nsample = train[(train.meter_reading < 5000)&(train.meter_reading > 0)] \nsample.sample(n=1000, random_state=1).hist(bins=50, figsize=(20,15)) ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# output a scatter matrix to see correlation between variables\n# meter reading is the target variable\n# the independent variables differ between each\nfrom pandas.plotting import scatter_matrix\nattributes = ['meter_reading', 'building_id', 'meter', 'site_id', 'primary_use']\nscatter_matrix(sample[attributes].sample(n=100, random_state=1), figsize=(12, 8))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"attributes = ['meter_reading','square_feet', 'year_built', 'floor_count']\nscatter_matrix(sample[attributes].sample(n=100, random_state=1), figsize=(12, 8))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"attributes = ['meter_reading', 'air_temperature', 'cloud_coverage', 'dew_temperature']\nscatter_matrix(sample[attributes].sample(n=100, random_state=1), figsize=(12, 8))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"attributes = ['meter_reading', 'sea_level_pressure', 'wind_direction', 'wind_speed']\nscatter_matrix(sample[attributes].sample(n=100, random_state=1), figsize=(12, 8))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The above illustrations shows that the square_feet column has right skew."},{"metadata":{"trusted":true},"cell_type":"code","source":"# the distribution of a sample from the square-feet column\ns = train['square_feet'].dropna(axis=0)\ns = s.sample(n=1000000, random_state=1)\nsns.distplot(a=s, kde=True)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can cube the data to get us closer a normal distribution (better than log(1+x)) and then replace the original values with the transformation."},{"metadata":{"trusted":true},"cell_type":"code","source":"# the distribution of the data behind the illustration cubed\ns = np.cbrt(s)\nsns.distplot(a=s, kde=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#transform column by cubing all values\ntrain.square_feet = np.cbrt(train.square_feet)\ntest.square_feet = np.cbrt(test.square_feet)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Feature Selection"},{"metadata":{},"cell_type":"markdown","source":"Sea-pressure lowers with altitdue. More energy is needed to heat air that is less dense. We divide air_temp by sea_level so that lower sea level pressures translates to a higher output. We square the results to remove the negative values."},{"metadata":{"trusted":true},"cell_type":"code","source":"train['sea_temp'] = np.sqrt(np.square(train['air_temperature'] / train.sea_level_pressure * 100))\ntest['sea_temp'] = np.sqrt(np.square(test['air_temperature'] / test.sea_level_pressure * 100))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Here we will try to approximate the air temperature that best regulates the building temperature without energy use in colder months. Sun exposure and heat produced naturally within the building (electronics, humans, etc. ) will lower the air temp needed to regulate a building. Eighteen degrees celcius is an comfortable indoor temperature. So, to get at a true zero temp we will shift the measurements down by 22 degrees then multiply these numbers by square feet so that any variance will be exagerrated. "},{"metadata":{"trusted":true},"cell_type":"code","source":"train['internal'] = np.square(train.air_temperature - 22)\ntest['internal'] = np.square(test.air_temperature - 22)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We divide the new variable, internal temperature, by the square feet of the building to approximate a ratio between building size and temperature. We square the results to remove negative values."},{"metadata":{"trusted":true},"cell_type":"code","source":"train['air_sq'] = np.square(train.internal / train.square_feet)\ntest['air_sq'] = np.square(test.internal / test.square_feet)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train.columns","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We will take a look at the correlation of these new columns to the target variable (meter_reading). "},{"metadata":{"trusted":true},"cell_type":"code","source":"# update the sample variable to include the new features\nsample = train[(train.meter_reading < 5000)&(train.meter_reading > 0)] \n\n\nattributes = ['meter_reading', 'sea_temp', 'air_sq']\nscatter_matrix(sample[attributes].sample(n=10000, random_state=1), figsize=(12, 8))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Drop Features"},{"metadata":{},"cell_type":"markdown","source":"The buildings represented in the dataset are located around the world (15 site_ids). This means that meterological and astronomical seasons of the northern and souther hemispheres differ by two seasons (summer = winter). Seasons, months, days, and years are inherently linked, so we will remove them. Hours and weekend column provide no additional correlation between readings. Other columns with a large number of missing values are also dropped. "},{"metadata":{"trusted":true},"cell_type":"code","source":"drop_train = ['floor_count', 'timestamp', 'year_built', 'cloud_coverage', 'precip_depth_1_hr']\ndrop_test = ['floor_count', 'timestamp', 'year_built', 'cloud_coverage', 'precip_depth_1_hr']\n\ntrain.drop(drop_train, axis=1, inplace=True)\ntest.drop(drop_test, axis=1, inplace=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(train.shape)\nprint(test.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Model"},{"metadata":{"trusted":true},"cell_type":"code","source":"#X = train.drop(columns=['meter_reading'])\n#y = train['meter_reading']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#from sklearn.model_selection import train_test_split, GridSearchCV\n#X_train, X_test, y_train, y_test = train_test_split(X,y, test_size=0.1, random_state=1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#categorical_features = [\"building_id\", \"site_id\", \"meter\", \"primary_use\", \"hour\", \"weekday\"]\n\n#lgb_train = lgb.Dataset(X_train, label=y_train, categorical_feature=categorical_features, free_raw_data=False)\n\n#lgb_test = lgb.Dataset(X_test, label=y_test, categorical_feature=categorical_features, free_raw_data=False)\n\nparams = {\n    \"objective\": \"regression\",\n    \"boosting\": \"gbdt\",\n    \"num_leaves\": 30,\n    \"learning_rate\": 0.1,\n    \"feature_fraction\": 0.85,\n    \"reg_lambda\": 2,\n    \"metric\": \"rmse\"\n}\n\n#model = lgb.train(params, train_set=lgb_train,  num_boost_round=1000, valid_sets=[lgb_train, lgb_test], verbose_eval=200, early_stopping_rounds=200)\n\n#predictions = model.predict(X_test, num_iteration=model.best_iteration)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#lgb_train_full = lgb.Dataset(X, label=y, categorical_feature=categorical_features, free_raw_data=False)\n\nparams = {\n    \"objective\": \"regression\",\n    \"boosting\": \"gbdt\",\n    \"num_leaves\": 40,\n    \"learning_rate\": 0.1,\n    \"feature_fraction\": 0.85,\n    \"reg_lambda\": 2,\n    \"metric\": \"rmse\"\n}\n\n#model2 = lgb.train(params, train_set=lgb_train_full,  num_boost_round=2000)\n\n#predictions = model2.predict(test, num_iteration=model.best_iteration)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Submission"},{"metadata":{"trusted":true},"cell_type":"code","source":"#submission = pd.DataFrame({'row_id':row_id, 'meter_reading':predictions})\n#submission.loc[submission.meter_reading < 0, 'meter_reading'] = 0\n#submission.to_csv('/Users/DataScience/energy/lgbm4.csv', index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#pd.read_csv('/Users/DataScience/energy/lgbm2.csv').head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#submission.head(10)","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"codemirror_mode":{"name":"ipython","version":3},"file_extension":".py","mimetype":"text/x-python","name":"python","nbconvert_exporter":"python","pygments_lexer":"ipython3","version":"3.7.3"}},"nbformat":4,"nbformat_minor":1}