Profiling, Benchmarking, and Simulation in R

“Programmers waste enormous amounts of time thinking about, or worrying about, the speed of noncritical parts of their programs, and these attempts at efficiency actually have a strong negative impact when debugging and maintenance are considered.” — Donald Knuth.

Optimising code to make it run faster is an iterative process:

  • Find the biggest bottleneck (the slowest part of your code).
  • Try to eliminate it (you may not succeed but that’s ok).
  • Repeat until your code is “fast enough.”

This sounds easy, but it’s not. (from Advanced R by Hadley Wickham)

Diagnostics and testing

First, there are many ideas on how to diagnose bottlenecks and improve performance. Reread what was just

  • where are bottlenecks: system.time() & proc.time()
  • what is broken: traceback() and debug()
  • does it do what I want?
    • informal unit testing of small pieces of code
    • using testthat package

Several issues of code efficiency come up:

  • cautions on using for and while loops
  • comparing floating point numbers: all.equal() and 1L
  • R profiling of code with Rprof() (before optimizing)

Example of system.time

surveys <- read.csv("../data/surveys.csv")
## Error in `file()`:
## ! cannot open the connection
forloop <- function(surveys) {
  n_species <- length(unique(surveys$species_id))
  means <- numeric(n_species)
  names(means) <- sort(unique(surveys$species_id))
  for(i in names(means))
    means[i] <- mean(surveys$weight[surveys$species_id == i], na.rm = TRUE)
  means
}
system.time(fmeans <- forloop(surveys))
## Error:
## ! object 'surveys' not found
library(dplyr)
dplyrmeans <- function(surveys) {
  surveys %>%
    group_by(species_id) %>%
    summarize(means = mean(weight, na.rm = TRUE)) %>%
    ungroup
}
system.time(dmeans <- dplyrmeans(surveys))
## Error:
## ! object 'surveys' not found

Alternatively, you can use the proc.time() function. See data example for timing using this approach.

It is a good idea to examine results from different approaches. Here are two ways to do this.

summary(fmeans - dmeans$means)
## Error:
## ! object 'fmeans' not found
all.equal(fmeans, dmeans$means)
## Error:
## ! object 'fmeans' not found

Example of traceback and debug

Suppose you have a little function that does not work quite right.

lousy <- function(x) {
  x <- as.character(x)
  y <- 1
  y <- sum(x, y)
  y
}
double_lousy <- function(x) {
  lousy(x)
}
double_lousy("a")
debug(lousy)

Do the following in the console. Step through lousy using return or n & return. You can examine parameters at each step, or try out the next step before running it.

double_lousy("a")

Unit tests

See Karl Broman’s Writing Tests. The idea is to write tests of units of code (unit tests), to check out code a step at a time. Then go further, and set up unit tests to be reused whenever you change your code.

Travis Continuous Integration works with packages on github to make sure code works as expected each time it is changed. For R, this is often used in conjunction with testthat. However, Travis/CI works for a wide range of software languages and platforms.


Simulation Studies in R

Learning Objectives:

  • Set up and run a statistical simulation (e.g. for a linear model).
  • Explain the value of simulation when analytical theory is unavailable.
  • Benchmark and compare computational speed of different simulation designs.

Simulation Lessons from Doug Bates (Stat 692)


This site uses Just the Docs, a documentation theme for Jekyll.