Enhanced Headroom Plot in R

September, 2008
by Adrian Cockcroft

About the Author
Adrian Cockcroft

Adrian is best known as the author of four books including Sun Performance and Tuning (2 editions); Resource Management; and Capacity Planning for Internet Services. In his 16 years at Sun he worked in technical sales and marketing, led creation of the BluePrints best practice publishing program, tested very complex integrated systems, was a leader of Sun's Six Sigma program and was the Chief Architect and Product Boss for Sun's High Performance Technical Computing business unit. In this time he gave many training classes and consulted with a wide range of customers, most notably as the on-site capacity planning consultant for the Salt Lake 2002 and Athens 2004 Olympic Games.

Joining eBay in 2004, he initially worked for Operations Architecture, investigating new platforms and providing guidance to the capacity planning groups at eBay and PayPal. As a founding member of eBay Research Labs in 2005, Adrian helped define the initial strategy for the Labs and an Innovation Forum. He researched operations related platforms and processes, lead research into advanced Skype plugin applications, contributed to development of the Skype4Java API and prototyped advanced wireless/mobile applications. During 2006 he published an IEEE paper on simulating large scale peer to peer networks, and a CMG paper on utilization measurement problems.

Adrian has consulted on architecture, scalability and performance for the Bebo.com social network, and is an advisory board member for Infovell and Holocosmos.

In 2007 Adrian joined Netflix as a Director of Web Engineering, directing a team responsible for research and development of scalable personalized web architectures.

Adrian filed two patents on capacity planning techniques while at Sun, and four patents related to peer to peer marketplaces while at eBay.

Adrian has a blog at http://perfcap.blogspot.com where he discusses capacity planning techniques, new computer technology, and how markets and innovation

interact. He is also a member of the Homebrew Mobile Phone Club, and several local classic car clubs.

For some reason I seem to find time to write code in R when I'm on an airplane. The last two trips I made resulted in significant enhancements and debugging of the code for my headroom plot. It started off simple but it now has a lot of bells and whistles, including color coding.

Main changes: the quantile used to remove outliers now only removes outliers that exceed the 95th percentile response time by default. It keeps all the throughput values unless you use qx=True.

In each of the throughput bins used to draw the histogram, the maximum response time for that bin is now calculated and displayed as a staircase line unless you set max=False.

The set of data is now split into ranges and color coded. The times series plot is coded so you can see the split, and the scatterplot shows how those points fall. I have been plotting weekly data at one minute intervals with split=7, which looks pretty good.

I read in some data that has been extracted from vxstat into a csv format at a known URL and plotted it three ways.

I plot the first 2880 data points, picking the read data rather than write, two days at one minute intervals.

> stime <- read.csv(url("http://somewhere/vxstat_response")) > sops <- read.csv(url("http://somewhere/vxstat_throughput")) > names(sops) [1] "DateTime" "vxstat_dg_operationsRead" [3] "vxstat_dg_operationsWrite" > chp(sops[1:2880,2],stime[1:2880,2]) > chp(sops[1:2880,2],stime[1:2880,2],q=1.0) > chp(sops[1:2880,2],stime[1:2880,2],q=1.0,splits=8)

Headroom Plot

Headroom Plot

Headroom Plot

Here is the code that generates the plot.

> chp <-function(throughput,response, q=0.95, qx=F, xl="Throughput",yl="Response", tl="Throughput Over Time",ml="Headroom Plot", fit=T, max=T, splits=0) { # remove zero throughput and response values nonzer <- (throughput != 0) & (response != 0) # array of true/false y <- response[nonzer] x <- throughput[nonzer] # remove outliers, keep response time points inside 95% by default if (q != 1.0) { quant <- (y < quantile(y,q)) # optionally trim throughput outliers as well if (qx) quant <- quant & (x < quantile(x, q)) x <- x[quant] y <- y[quant] } # make histograms and record end points for scaling xhist <- hist(x,plot=FALSE) yhist <- hist(y,plot=FALSE) xbf <- xhist$breaks[1] # first ybf <- yhist$breaks[1] # first xbl <- xhist$breaks[length(xhist$breaks)] # last ybl <- yhist$breaks[length(yhist$breaks)] # last xcl <- length(xhist$counts) # count length ycl <- length(yhist$counts) # count length xrange <- c(0.0,xbl) yrange <- c(0.0,ybl) xlen <- length(x) # make a multi-region layout nf <- layout(matrix(c(1,3,4,2),2,2,byrow=TRUE), c(3,1), c(1,3), TRUE) layout.show(nf) # set plot margins for throughput histogram and plot it par(mar=c(0,4,3,0)) barplot(xhist$counts, axes=FALSE, xlim=c(xcl*0.00-xbf/((xbl-xbf)/(xcl-0.5)),xcl*1.00), ylim=c(0, max(xhist$counts)), space=0, main=ml) # set plot margins for response histogram and plot it sideways par(mar=c(5,0,0,1)) barplot(yhist$counts, axes=FALSE, xlim=c(0,max(yhist$counts)), ylim=c(ycl*0.00-ybf/((ybl-ybf)/(ycl-0.5)),ycl*1.00), space=0, horiz=TRUE) # set plot margins for time series plot par(mar=c(2.5,1.7,3,1)) plot(x, main=tl, cex.axis=0.8, cex.main=0.8, type="S") if (splits > 0) { step <- xlen/splits for(n in 0:(splits-1)) { lines((1+n*step):min((n+1)*step,xlen), x[(1+n*step):min((n+1)*step,xlen)], col=4+n) } } # set plot margins for main plot area par(mar=c(5,4,0,0)) plot(x, y, xlim=xrange, ylim=yrange, xlab=xl, ylab=yl, pch=20) if (max) { # max curve b <- xhist$breaks i <- b[2] - b[1] # interval maxl <- list(y[b[1] < x & x <= (b[1]+i)]) for(n in b[c(-1,-length(b))]) maxl <- c(maxl,list(y[n < x & x <= (n+i)])) #print(maxl) maxv <- unlist(lapply(maxl,max)) # apply max function to elements of list #print(maxv) #lines(xhist$mids,maxv,col=2) # join the dots #staircase plot showing the range for each max response lines(rep(b,1,each=2)[2:(2*length(maxv)+1)],rep(maxv,1,each=2),col=3) } if (fit) { # fit curve, weighted to predict high throughput # create persistent chpfit object using <<- chpfit <- glm(y ~ x, inverse.gaussian, weights=as.numeric(x)) # add fitted values to plot, sorted by throughput lines(x[order(x)],chpfit$fitted.values[order(x)],col=2) } if (splits > 0) { step <- xlen/splits for(n in 0:(splits-1)) { Sys.sleep(1) points(x[(1+n*step):min((n+1)*step,xlen)],y[(1+n*step):min((n+1)*step,xlen)], xlim=xrange, ylim=yrange, col=4+n) } } }