{"cells":[{"metadata":{},"cell_type":"markdown","source":"## ASHRAE Analysis\n\nThe idea for this notebook is to explore a few important questions about this competition. I'll be updating it for the next few weeks.\n\n\n### Our target\n\nThis competition is a little different: we are predicting energy consumption backwards (for the past). But why are we doing that? \n\nThe key idea here is that we have some consumption history for each building, but in a given point we invested in new equipament to reduce the energy cost for the building. After this invesment, we need to evaluate how much the building is salving on energy, so a naive approach would be to compare the current consumption (after the improvements) with past consumption - using some mean for instance. This would be ok, but we can do better than that. So in this competition we use a model to predict the consumption (without new equipament) based on the weather and other important conditions to make this comparison more accurate.\n\n![](https://www.mdpi.com/make/make-01-00056/article_deploy/html/images/make-01-00056-g001-550.jpg)\n<font size=\"1\"><center>Image from [this notebook](https://www.kaggle.com/hmendonca/starter-eda-and-feature-selection-ashrae3)</center></font>"},{"metadata":{},"cell_type":"markdown","source":"## Overview\n\n* Training data has 16 columns and 20.2 MM rows from  2016-01-01 to 2016-12-31\n* Test data has 41.7 MM rows from 2017-01-01 to 2018-05-09\n* There are 1,449 different buildings in 16 locations (site_id)\n* Buildings and sites are the same for train and test"},{"metadata":{"_uuid":"051d70d956493feee0c6d64651c6a088724dca2a","_execution_state":"idle","trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"library(tidyverse)\nlibrary(repr)\nlibrary(lubridate)\noptions(tibble.width = Inf)\noptions(warn=-1)\n\n# theme for ggplot\ntheme_set(theme_minimal() +\n            theme(axis.title.x = element_text(size = 12, hjust = 1),\n                  axis.title.y = element_text(size = 12),\n                  axis.text.x = element_text(size = 11),\n                  axis.text.y = element_text(size = 11),\n                  panel.grid.major = element_line(linetype = 2),\n                  panel.grid.minor = element_line(linetype = 2),\n                  plot.title = element_text(size = 18, colour = \"grey25\", face = \"bold\")))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dpath <- \"../input/ashrae-energy-prediction/\"\nprint(list.files(path = dpath))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# read the data\ntrain <- read_csv(paste(dpath, \"train.csv\", sep=\"\"), col_types=cols())\nbuilding <- read_csv(paste(dpath, \"building_metadata.csv\", sep=\"\"), col_types=cols())\nweather_train <- read_csv(paste(dpath, \"weather_train.csv\", sep=\"\"), col_types=cols())\n\n# join train data\ntrain <- train %>%\n  left_join(building, by = \"building_id\") %>%\n  left_join(weather_train, by = c(\"site_id\", \"timestamp\"))\nrm(weather_train)\nprint(train, n=5)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"summary(train)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Meter and Meter Reading"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"options(repr.plot.width=12, repr.plot.height=5)\nggplot(train) + geom_histogram(aes(x=meter_reading), bins=40) +\nfacet_wrap(~ factor(meter), nrow=2, scales=\"free\") + ggtitle(\"Reading histogram\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* There are outliers in our target variable for all four meters"},{"metadata":{},"cell_type":"markdown","source":"There are four meter types:\n0. electricity\n1. chilledwater\n2. steam\n3. hotwater"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"ggplot(train) + geom_bar(aes(x=factor(site_id), fill=factor(meter)))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* Most data are for eletricity consumption\n* All sites have an eletricity meter\n* Sites 3, 4, 5, 8 and 12 have only eletricity\n* Hotwater is a small fraction of the data"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"train$date <- as.Date(format(train$timestamp, \"%Y-%m-%d\"))\ntotal_reading <- train %>% group_by(date, meter) %>% summarize(meter_reading=median(meter_reading, na.rm=TRUE))\nggplot(total_reading) +\ngeom_smooth(aes(x=date, y=meter_reading, color=factor(meter)), method='gam', formula = y ~ s(x, bs = \"cs\")) +\nggtitle(\"Median daily reading\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* We have quite different readings depending on the meter (scale problem)\n* Using one model for each meter might be a good idea\n* Seasonality: more eletricity/chilledwater on summer and steam/hotwater on winter\n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"options(repr.plot.width=12, repr.plot.height=8)\ntotal_reading <- train %>% group_by(site_id, date, meter) %>% summarize(meter_reading=median(meter_reading, na.rm=TRUE))\nggplot(total_reading) +\ngeom_smooth(aes(x=date, y=meter_reading, color=factor(site_id)), method = 'gam', formula = y ~ s(x, bs = \"cs\")) +\nfacet_wrap(~ factor(meter), nrow=2, scales=\"free\") + ggtitle(\"Median daily reading for each site and meter\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* Some sites have much higher readings, which could be a consequence of the weather or bulding size (area/floors)\n* It seems that the locations are at the northern hemisphere"},{"metadata":{},"cell_type":"markdown","source":"## Weather\n\nLet's check the correlation between temperature and meter reading:"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"options(repr.plot.width=14, repr.plot.height=8)\nair_temp <- train %>% group_by(timestamp, site_id, meter) %>% summarize(air_temperature=mean(air_temperature, na.rm=TRUE), meter_reading=median(meter_reading))\nggplot(filter(air_temp, meter_reading < 10000)) +\ngeom_point(aes(x=air_temperature, y=meter_reading, color=factor(site_id))) +\nfacet_wrap(~ factor(meter), nrow=2, scales=\"free\") +\nggtitle(\"Median reading for air temperature\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"dew_temp <- train %>% group_by(timestamp, site_id, meter) %>% summarize(dew_temperature=mean(dew_temperature, na.rm=TRUE), meter_reading=median(meter_reading))\nggplot(filter(dew_temp, meter_reading < 10000)) +\ngeom_smooth(aes(x=dew_temperature, y=meter_reading, color=factor(site_id)), method = 'gam', formula = y ~ s(x, bs = \"cs\")) +\nfacet_wrap(~ factor(meter), nrow=2, scales=\"free\") +\nggtitle(\"Median reading for dew temperature\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* Temperature is more relevant to some sites than others\n* Low impact for eletricity\n* Air and dew temperature seems quite similar"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"speed <- train %>% group_by(timestamp, site_id, meter) %>% summarize(wind_speed=mean(wind_speed, na.rm=TRUE), meter_reading=median(meter_reading))\nggplot(filter(speed, meter_reading < 10000)) +\ngeom_smooth(aes(x=wind_speed, y=meter_reading), method = 'gam', formula = y ~ s(x, bs = \"cs\")) +\nfacet_wrap(~ factor(meter), nrow=2, scales=\"free\") +\nggtitle(\"Median reading for wind speed\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Time analysis"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"train$hour <- hour(train$timestamp)\noptions(repr.plot.width=14, repr.plot.height=12)\nhour_data <- train %>% group_by(hour, site_id) %>% summarize(air_temperature=median(air_temperature, na.rm=TRUE), dew_temperature=median(dew_temperature, na.rm=TRUE))\nggplot(hour_data) +\ngeom_smooth(aes(x=hour, y=air_temperature), method = 'gam', formula = y ~ s(x, bs = \"cs\")) +\nfacet_wrap(~ factor(site_id), scales=\"free\") +\nggtitle(\"Median air temperature for each hour of the day and site\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* Air temperature and hour doesn't match"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"train$month <- month(train$timestamp)\nmonth_data <- train %>% group_by(month, site_id) %>% summarize(air_temperature=median(air_temperature, na.rm=TRUE), dew_temperature=median(dew_temperature, na.rm=TRUE))\nggplot(month_data) +\ngeom_smooth(aes(x=month, y=air_temperature), method = 'gam', formula = y ~ s(x, bs = \"cs\")) +\nfacet_wrap(~ factor(site_id), scales=\"free\") +\nggtitle(\"Median air temperature for month and site\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* It seems that all sites are in the Northern Hemisphere"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"month_data <- train %>% group_by(month, site_id, meter) %>% summarize(meter_reading=median(meter_reading, na.rm=TRUE))\nggplot(month_data) +\ngeom_smooth(aes(x=month, y=meter_reading, color=factor(meter)), method = 'gam', formula = y ~ s(x, bs = \"cs\")) +\nfacet_wrap(~ factor(site_id), scales=\"free\") +\nggtitle(\"Meter reading for month and site\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Building"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"options(repr.plot.width=14, repr.plot.height=8)\nsqft <- train %>% group_by(building_id, year_built, meter) %>%\nsummarize(square_feet=mean(square_feet, na.rm=TRUE), meter_reading=median(meter_reading))\n\nggplot(sqft) +\ngeom_point(aes(x=square_feet, y=meter_reading, color=year_built)) +\nfacet_wrap(~ factor(meter), nrow=2, scales=\"free\") +\nggtitle(\"Median reading by area\")","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"options(repr.plot.width=12, repr.plot.height=10)\nggplot(building) + geom_bar(aes(y=factor(primary_use), fill=factor(site_id))) +\nggtitle(\"Number of buildings by Primary use\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"* Most buildings are for education\n* Very few data for some categories, which might be a problem since they are very different (eg. manufacturing vs office)\n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"options(repr.plot.width=14, repr.plot.height=12)\nggplot(building) + geom_point(aes(x=square_feet, y=floor_count, color=year_built)) +\nfacet_wrap(~ factor(primary_use)) + ggtitle(\"Area and Floor count by primary use\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Working in progress.."}],"metadata":{"kernelspec":{"display_name":"R","language":"R","name":"ir"},"language_info":{"mimetype":"text/x-r-source","name":"R","pygments_lexer":"r","version":"3.4.2","file_extension":".r","codemirror_mode":"r"}},"nbformat":4,"nbformat_minor":1}