Author Archives: Charles_Mohan

Visualizing a Power Buying Strategy in R

In the retail energy space, customers often have the option to fix their future energy costs at current market prices, or fix their energy costs over time at varying prices. Sometimes it can be difficult to communicate / visualize the difference between these two options. This post will show one to do this and also provide an excuse to play with the awesome animation package in R!

Data: Here Code: Here

The chart below is an example of a completely naive buying strategy – buy 1/12 of your power needs on the 15th of each month.  The data below is the 12 month average ATC (around the clock) PJM Western Hub price for Jul13 – Jun14 and Jul14 – Jun15.  This would result in the 2013-2014 price being fixed by July 2013, and 2014-2015 price being fixed by July 2014. This would produce a fixed power price per megawatt by the start of the fiscal year, and help avoid budget uncertainty due to fluctuations in energy prices.

buyChart

The large dots represent buys and the horizontal line is the final average price per megawatt for each year.  Now, this chart shows a lot of information, but to really drive home how this strategy will work over time I’ve used the animation package in R to ‘grow’ the chart. The final result is the video below:

 

As you can see, there is additional information conveyed when you see the above chart evolve over time that you don’t see in the final slide. Particularly, you get to see the horizontal lines move over the course of the year – visualizing the ‘path’ to the final fixed price or, depending on the strategy, the amount of volatility that should be expected.

The main package used to create the above video is the animation package in R. The resulting .mp4 video can be uploaded to YouTube as I have done, or if needed for a presentation, can be easily embedded into a slide deck.

library(dplyr)
library(lubridate)
library(ggplot2)
library(animation)
library(tidyr)
library(RCurl)

First we need to get a set of dates on which we will buy, for this example we are picking the closet weekday to the 15th of each month. This will obviously change depending on the strategy being employed.

buys<-c('7/13/2012','8/15/2012','9/14/2012',
'10/15/2012','11/15/2012','12/14/2012',
'1/15/2013','2/15/2013','3/15/2013',
'4/15/2013','5/15/2013','6/14/2013',
'7/15/2013','8/15/2013','9/13/2013',
'10/15/2013','11/15/2013','12/13/2013',
'1/15/2014','2/14/2014','3/14/2014',
'4/15/2014','5/15/2014','6/13/2014'
)

Download the data from GitHub and do some data munging.

url <- "https://raw.githubusercontent.com/cam333/blog_posts/master/blogCode_animationChart/13_14.csv"
x <- getURL(url)

rawdata <- read.csv(text = x) %>%
select(Date,July_13_June_14_fw,July_14_June_15_fw) %>%
mutate(Date = mdy(Date)) %>%
filter(Date > ymd('2012-06-29')) %>%
mutate(buy = as.character(ifelse(Date %in% mdy(buys),1,0)),
July_13_June_14_fw = July_13_June_14_fw,
July_14_June_15_fw = July_14_June_15_fw) %>%
gather(variable,value,-Date,-buy) %>%
na.omit()

Next we need to write a function that creates our chart. Notice how the labels, point formatting, and hline geom are all dependent on attributes of the data being charted. This is important since the animation will be created by generating many charts, each with a slightly larger data set – allowing the labels to adjust as we move through time.

chart<-function(cutoff){

a <-ggplot(filter(rawdata,Date <= cutoff),aes(x=Date, y=value, group = variable)) +
geom_line() +
geom_hline(yintercept = filter(rawdata,Date <= cutoff &
buy == 1 &
variable == 'July_13_June_14_fw') %>%
summarize(mean = mean(value, na.rm=TRUE)) %>% use_series(.,mean),
colour = "red", geom_text(aes(label = "2013 Avg."))
) +
geom_hline(yintercept = filter(rawdata,Date <= cutoff &
buy == 1 &
variable == 'July_14_June_15_fw') %>%
summarize(mean = mean(value, na.rm=TRUE)) %>% use_series(.,mean),
colour = "blue", geom_text(aes(label = "2013 Avg."))
) +
geom_point(aes(size=buy,alpha = buy, color = variable, fill = buy)) +
scale_colour_manual(values = c("red","blue")) +
ylim(35,50) +
xlim( ymd('2012-06-29'),ymd('2014-06-30')) +
annotate("text", x =ymd('2013-12-01'), y = 50,
label = paste("Average Price 2013-2014: $", as.character(round(filter(rawdata,Date <= cutoff &
buy == 1 &
variable == 'July_13_June_14_fw') %>%
summarize(mean(value, na.rm=TRUE)),
digits=2)),sep =" "),colour = 'red',hjust =0) +
annotate("text", x =ymd('2013-12-01'), y = 49,
label = paste("Average Price 2014-2015: $", as.character(round(filter(rawdata,Date <= cutoff &
buy == 1 &
variable == 'July_14_June_15_fw') %>%
summarize(mean(value, na.rm=TRUE)),
digits=2)),sep =" "), colour = 'blue', fontface = 'bold',hjust =0)+
annotate("text", x =ymd('2013-12-01'), y = 48,
label = paste("Date: ",filter(rawdata,Date <= cutoff) %>% summarize(date = last(Date))%>%
extract2(.,1),sep=" " ),hjust =0, fontface = 'bold') +

theme_bw(base_size = 15) +
theme(legend.position="none") +
ggtitle("Simulation of Monthly Buying") +
ylab("$/MWh") +

geom_point(data = data.frame(x=ymd('2012-07-29 UTC'),y=50,variable = 'July_13_June_14_fw'),aes(x,y),color = 'red',size=7) +
annotate("text", x =ymd('2012-11-01 UTC'), y = 50,label = "Buy Points 2013-2014") +

geom_point(data = data.frame(x=ymd('2012-07-29 UTC'),y=48,variable = 'July_14_June_15_fw'),aes(x,y),color = 'blue',size=7) +
annotate("text", x =ymd('2012-11-01 UTC'), y = 48,label = "Buy Points 2014-2015")

print(a)

}

Now we need a function that loops over a sequence of dates and produces a chart each time. This is the function we will call with the ‘saveVideo’ function in the animation package.

animate <- function() {
lapply(seq(ymd('2012-06-29 UTC'),ymd('2014-06-30'),by = '2 day'), function(i) {
chart(ymd(i))
})
}

Lastly we need to run the animate function and save the video. This can take some time depending on how complicated your chart is and how many times you are updating it. We will be producing 366 charts to create our video.

saveVideo(animate(),interval=.05, video.name="C:\\Users\\cmohan\\Desktop\\buypoint2.mp4",
ani.width = 800,
ani.height = 372)

You should now have a video clip on your desktop.

Creating Hourly Electricity Price Scalars in R

Forward prices in the power markets are typically quoted on a monthly basis. However, it is often helpful to construct an hourly price series (whose average equals the average for the month) for use in risk management, hedging strategy decisions, and retail pricing. This post will show one way this can be done.

A simple way to construct an hourly price series from a monthly price, is to use a price scalar. This is nothing more than the historical ratio associated with a particular month and hour. For example, it answers the question: “Historically, how does hour 15 in January relate to January’s monthly average?”  Note: 1) This is not a price forecast in the traditional sense – think of it as a way to estimate typical hourly deviations from the average., and 2) this is one of the most basic ways to make this estimate.     

The below steps will show you how to answer the above question relying heavily on the dplyr and lubridate packages. GitHub code here.

library(dplyr)
library(lubridate)
library(ggplot2)
library(tidyr)

The Data:

First, we need to get a hold of historical price data.  I’m using the PJM total LMP Western Hub price accessible from their Data Miner service.  The options I’ve selected are shown in the screen shot below:

dataminer_options

Once you have downloaded the data as a .csv file, it is time to do some data munging. The code below takes the data from wide format to long, and fixes the date to include hours.

raw_data <- read.csv('C:\\Users\\cmohan\\Desktop\\DataMiner-Export_2016-01-12-westhub.csv') %>%
select(-VERSION,-ZONE,-PNODEID,-PNODETYPE,-PRICINGTYPE,-H2DST) %>%
mutate(PUBLISHDATE = mdy(PUBLISHDATE)) %>%
gather(hour,value,-PUBLISHDATE,-PNODENAME) %>%
mutate(variable = as.numeric(substring(hour,2)) - 1,
Date = mdy_hms(paste(month(PUBLISHDATE),day(PUBLISHDATE),year(PUBLISHDATE),hour,0,0))) %>%
select(Date,value) %>%
arrange(Date) %>%
na.omit()

Now that the data is in a form we can work with, it’s time to think about what we are trying to do. The point of an hourly scalar is to take the monthly forward price and extrapolate it to hourly granularity. To do this we will need a scalar that increases or decreases the monthly forward price based on the historical relationship between the monthly average price and each particular hour.

I am focusing on three subsets or types of hourly power scalars: On Peak, Off Peak, and Weekends.

powerhours

With this type of grouping, each hour will have two different scalars depending on whether the day is in the workweek or on the weekend (or holiday).

Calculating the Scalar:

So, the first thing we need to do is add a ‘type’ column to the historical price data to identify each hour as either On Peak, Off Peak, or Weekend. The below code will create this type column.

 filtered_data <- raw_data %>%
mutate(value = as.numeric(value),
type = ".",
type = ifelse(wday(Date) %in% seq(2,6) & hour(Date) %in% seq(7,22),'OnPeak','OffPeak'),
type = ifelse(wday(Date) %in% c(1,7), 'WeekEnd',type)
)

types <- c('OnPeak','OffPeak','WeekEnd')

 

Now we have each hour in the historical data set labeled as OnPeak, OffPeak, or Weekend (For simplicity, holidays are being ignored).  The next step is to  1) split the dataframe by year, 2) calculate the average of each month, and 3) divide each hour by the monthly average. These three steps need to be looped over for each month and type. (there’s probably a way to do this all with apply functions, suggestions are more than welcomed)

 out <- list()
# for each month
for(m in 1:12){
# for each type
for (t in 1:3 ){
fdata <- filter(filtered_data,type == types[t],month(Date) == m)

by_year <- split(fdata, year(fdata$Date))
pct_of_mean <- lapply(by_year, function(x) mutate(x,value = value/mean(value)))
s <- do.call(rbind, pct_of_mean)

name <- paste(types[t],':', month(m, label = TRUE),sep="")


out[[name]] <- s
}

}

all <- do.call(rbind,out)

 

The ‘all’ data.frame contains all the scalars by hour, all we need to do is summarize this by month, hour, and type. The code below does just that.

 

 scalars <- arrange(all) %>%
mutate(month = month(Date),hour = hour(Date))%>%
group_by(month,hour,type) %>%
summarize(scalar = mean(value)) %>%
ungroup()

 

Now we can chart the results.Rplot

 ggplot(mutate(scalars, month = as.factor(month), hour = as.factor(hour)), aes(x = hour, y = scalar,group = month)) +
geom_line(aes(colour = month)) +
scale_y_continuous(breaks = seq(-2,2,by = .25),'Scalar') +
facet_grid(type~.) +
theme_bw()

Using the scalar:

Now that we have scalars, let’s use them. Say we want to create an hourly shaped price for August 2016.  The first thing to do is create dates we want to forecast.

 dates <- data.frame(date = seq(ymd_h('2016-08-01-0'),ymd_h('2016-08-10-23'), by = 'hours')) %>%
mutate(month = month(date),
hour = hour(date),
type = ".",
type = ifelse(wday(date) %in% seq(2,6) & hour(date) %in% seq(7,22),'OnPeak','OffPeak'),
type = ifelse(wday(date) %in% c(1,7), 'WeekEnd',type))

And we need a forward price for West Hub in August.

  forwardPrice <- data.frame(month = 8, type = c('OnPeak','OffPeak','WeekEnd'),price = c(52,35,32))

With a few left joins we can easily get an hourly representation of the monthly price for the first 10 days in August 2016.

 

 fcast <- Reduce(left_join,list(dates,scalars, forwardPrice)) %>%
mutate(fcast_price = scalar * price)

ggplot(fcast, aes(x= date, y = fcast_price)) +
geom_line() +
theme_bw()

fcast_plot