I have hourly temperature and power consumption data of several days of a month. The pattern is almost similar across days like this: enter image description here

Using this data I want to predict the usage of a coming day. I have features : 1) hour of the day 2) temperature; and the response variable, power. Looking at the data, I believe I should fit three separate models and not a single model

  1. First model for data from Midnight to 10 AM as usage remains almost constant during this time, and temperature does not vary too much
  2. Second model for data from 11AM till 6 PM. This portion follows a sharp increase and then almost constant usage
  3. Third model for data from 7 PM till midnight. This portion shows constant decrease in power

To follow this intuition, I used three models accordingly and later combined predictions from these to output sequence of 24 numbers for coming day. Formula for each of these models is: lm(power ~ time_hour + temperature, data = xxx), but each model is trained with the data corresponding to specific time duration of the day.

Instead of dividing the data manually and using three separate models, is there any other existing technique which will take care of our intuition and does not need manual division of data or creation of separate models.

During my search, I found that I can use GAM (Generalized Additive Models) and I came up with following formula

lm(power ~ ns(time_hour, knots =(9, 18)) + temperture, data = xxx)

Using above formula I think that I am placing knots at 9 AM and at 6 PM. Right? I don't know how should I enforce knots in temperature feature exactly at these specific times, so that knots of temperature and time_hour will sync.

The above plots were plotted using following data:

dframe <- structure(list(time_hour = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 
12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24), temperature = c(22.5, 
24, 23.5, 20.5, 22.5, 22.5, 19.5, 23.5, 23, 20.5, 26.5, 28.5, 
30, 32, 33.5, 33, 30.5, 30, 29.5, 29, 28, 27, 28, 28.5), power = c(97.04319, 
95.7225, 88.59191, 88.34882, 90.17179, 88.82062, 87.73833, 89.36342, 
85.31775, 91.1292, 116.79035, 149.58614, 172.32438, 171.27931, 
159.53858, 162.03544, 170.78468, 164.0275, 155.86717, 135.77197, 
133.01235, 116.29253, 100.87483, 97.84942)), .Names = c("time_hour", 
"temperature", "power"), row.names = c(NA, -24L), class = "data.frame")

and minimal code used is:

par(mfrow= c(1,2))
plot(dframe$time_hour ,dframe$power,type="l",xlab = "Hour of day", ylab = "power" )
plot(dframe$time_hour ,dframe$temperature,type="l",xlab = "Hour of day", ylab = "Temperature" )
  • $\begingroup$ if the answer that I posted works for you it would be appreciated if you accept the answer. $\endgroup$
    – Stereo
    Oct 23 '16 at 15:32

This is a typical time series problem. The first step is to make sure your time series is stationary, see here for a good explanation why this is necessary.

So we first do,

dframe$power_diff <- c(NA, diff(dframe$power))
plot(dframe$power_diff, type = "l")



Next we want to consider an ARIMAX model. As part of this exercise we want to know whether there is a dependency in the time series, a so called auto regressive term. We can research this using,


which gives,


In this plot we see that there is a positive auto regressive term at lag 1 at a significance level of 95%, we conclude this based on the observation that the second line, which is lag 1, is above the blue dotted line.

Next we want to fit a ARIMAX model.

# Convert data to ts object
power_data <- ts(data = dframe$power_diff[-1], frequency = 8766, start = c(2016, 1, 1))
temp_data <- ts(data = dframe$temp_diff[-1], frequency = 8766, start = c(2016, 1, 1))

# Build the model.
m <- Arima(power_data, order = c(1, 0, 0), xreg = temp_data)

# Let's see how our fit is.
lines(fitted(m), col = "blue")

This gives this result,


You can play with the ARIMA(p, i, q) terms. We should leave $p=1$ since we found that the AR term is significant at lag 1. If we add a MA(1) we find a model that has a MAPE=260.46 versus a MAPE=272.62 for our initial model. So you can create the best model using,

m2 <- Arima(power_data, order = c(1, 1, 0), xreg = temp_data)
lines(fitted(m2), col = 'red')

Will lead to, arimas

Of course we can automate the above step with auto.arima(). This function will search for the best arima model,

m3 <- auto.arima(power_data, xreg = temp_data)
lines(fitted(m3), col = 'green')

which will give you, auto_arima

You see that each model has a poor fit in the spike. To improve this you could introduce a dummy variable,

xdata <- ts(data = as.matrix(cbind(temp_data, dframe$power_diff[-1] > 20)), frequency = 8766, start =c(2016,1,1))
m_dummy <- auto.arima(power_data, xreg = xdata)
lines(fitted(m_dummy), col = 'orange')

which gives us the following result, dummy

We see that this model gives us a much better result. So if you can relate dummy to a point in time you could create a model that has a good fit.

According Occam's Razor you should try to keep your models simple. So I would suggest to build a single model for all your time slots instead of three separate models as you suggested in your question.

You can turn this back to the original time series with cumsum. Hope this helps.

  • $\begingroup$ I have some experience with SARIMA. My issue is that I have time series data which contains seasonality at day level and week level. When I use SARIMA it takes lots of time to get the results. Becuase of the time issue, I am trying to avoid SARIMA. $\endgroup$ Oct 24 '16 at 13:54
  • $\begingroup$ What takes a lot of time? If finding the right model is time consuming you may give seasonal a try. I would recommend to add this constraint to your question. Or vote for my question and post a new question that has your more extended question. $\endgroup$
    – Stereo
    Oct 24 '16 at 14:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.