Code for the NYC R Meetup on debugging R

Notes.txt

If you see the bug, keep your damn mouth shut!

Example #1
Show data.
I've written pointless code. 
traceback() shows where the problem lies, but not when or why
Explain "call stack" made of "frames", noting that many things that don't 
syntactically look like function calls in R actually are.
Error doesn't seem to happen all the time, so adding browser() calls 
might not be very helpful.
Use options(error=recover) to invoke recover, which is basically
a wrapper around browser, only when the error happens.
Don't know why frames are in opposite order in recover vs traceback!
Unlike browser, "c" or "n" or Enter drop you back to the stack trace.
Look at frame 1 first. Use ls() to see variables. head(df), head(df2),
tail(df2).
Lookat frame 2 now. See that i = 60. Use body() to see the body of the current
function. df[i, "wacky.sum"] looks ok. df[i, "rx"] looks ok. df[i, 'axyzoz.XS']
gives NULL! That's odd. Look at whole row with df[i, ]. Row name is wrong?
unique(df$rx), names(df). Yes.
Quit out, fire the intern, edit that CSV file, rerun, no problemo.
Turn off with options(error=NULL) 

Example #2
head(cats.test), head(rxs.test)
Walk through code. Merge rx effectiveness into pets df. Create a vector based
on the loss vector. Copy the df, add the new column, and return. (Yes, silly in
this case!)
Try it: uh oh
traceback() says the problem is in the assignment to exp.profit field, toward the end.
debug(do.tricky.thing)
debug is like browser, except you can walk through statements in the current block
with "n" or Enter. "c" drops you out of the debugger and return to normal execution.
"Q" quits entirely.
sim.result numbers look reasonable, but there are 7, and pets has only 6 rows!
merged.test has 7 rows two! baby is duplicated! something must be wrong with merge.
table(meds$rx) shows a dupe. That'll do it.
Assert that there are no dupes! That puts the error in the right place.
stopifnot(all(table(meds$rx) == 1))
(or, if (...) stop ("error msg"))
Confirm that that is the case.
Solve the problem by going back to the data source and correcting it.


Example #3
Show two data frames.
Walk quickly through code. Run it.
There's no error, do the numbers make sense?
Well, the reported numbers should be roughly the mean of the Ns for each medicine.
The one for Prlylife is too small! WTF?
If we had temporary variables in scope, we could print those, but everything's buried in a 
ddply in a function!
We might want to debug within the anonymous ddply function, but we don't want to 
edit the code source, and we don't want to have to walk through what might be thousands
of calls to that function just to get to the error!
Trick is called "conditional breakpoints". They're not that easy to set in R, sadly,
although I suspect they're easier with Revo R's GUI.
To use conditional breakpoints in R, use trace() function. Let's walk through it...
?trace
Most useful parameters are what, tracer, and at.
The tracer parameter is an expression. By default, it's a function that prints the call.
So, trace(sample) shows each call to the sample function. Demo.
You can also use other expressions, like the browser function. trace(sample, browser).
Demo. 
browser shows the original code
body() shows that the call to browser was added into the version of the function in memory!
Woah, super useful!
But this function's not immediately very useful to us, at least not right now.
What we really want is to invoke browser from within the anonymous function that's the
argument to ddply! But it has no name, so we can't! The at argument to trace is the 
solution. It's difficult to use, though. 
untrace(sample)
as.list(body(do.trickier.thing))
as.list(body(do.trickier.thing)[[3]])
as.list(body(do.trickier.thing)[[3]][[3]][[4]])
as.list(body(do.trickier.thing))
trace(do.trickier.thing, browser, at=list(c(3,3,4,3)))
do.trickier.thing(rxs.test, exps.test)
body()
Now in the right place, but we don't actually want to stop now, since this one seems
to work! Make this conditional:
trace(do.trickier.thing, quote(if (df$rx[[1]]=='prlylyfe') browser()), at=list(c(3,3,4,3)))
body() shows the new conditional code
df shows the data
Now let's try evaluating something, and finally fixing this bug.
Confirm call to data.frame gives answer about half of what's expected, which is a bit more than
16. 
Debug from the inside out, start with sample(df$Ns, 10, replace=TRUE). Wait, what? This should
give me 10 16s, but it doesn't! 
?sample
This is what you get when people try to be clever. Software library design is no time to be 
clever!
Many solutions. Might just use if(length(df$Ns) == 1) rep() or something...

piratical.R

options(stringsAsFactors=FALSE)
library(plyr)

# Example #1: Complex processing, don't want to modify code. traceback(),
# error=recover().

all.dat <- read.csv('piratical_dat.csv')
dat.test <- head(all.dat)

process.complexiously <- function(df, arg1 = .12, arg2 = function(x, q=arg1) ifelse(x<q, x, q+2), 
					    stuff=0) {
	do <- sum(stuff)

	other <- rep(arg1, 44)
	do <- sum(other, stuff)
	
	df2 <- ddply(df, .(rx), transform, wacky.sum=dzyhapi + qwexinox + 7)

	more.wackiness(df2)
}

more.wackiness <- function(df) {
	for (i in 1:(nrow(df) * .95)) 
	{
		df[i, 'wacky.sum'] <- df[i, 'wacky.sum'] + sum(df[, df[i, 'rx']])
	}
	df
}

process.complexiously(dat.test)

process.complexiously(all.dat)


# Example #2: Error doesn't tell you where the problem is. Debug and walk through. Asssert
# invariants.

cats.test <- data.frame(cat = c('fluffy', 'spot', 'tiger', 'baby', 'garth', 'felix'),
				height = c(15, 9, 22, 11, 12, 9),
				weight = c(6, 5, 14, 9, 4, 5),
				rx = c('dzyhapi', 'qqrvibe AB', 'xyxyslim', 'ytrwow', 'prlylyfe', 'dzyhapi'))


rxs.test <- data.frame(rx = c('dzyhapi', 'qqrvibe AB', 'aarvbyol', 'prlylyfe', 'ytrwow', 'zorncat',
					'qwexinox', 'ytrwow', 'xyxyslim', 'axyzoz XL', 'poyqwr'),
				loss = c(.04, .13, .02, .04, .06, .10, .09, .04, .01, .05, .05))

do.tricky.thing <- function(pets, meds, sd=.03)
{
  merged.test <- merge(pets, meds, by='rx')

	# complicated Monte Carlo simulation
	sim.result <- with(merged.test, loss + rnorm(length(loss), mean=0, sd=sd))

	# put it back together and return
	ret <- pets
	ret$exp.profit <- sim.result
	ret
}

res.test <- do.tricky.thing(cats.test, rxs.test)



# Example #3: Debugging when there's no error. Conditional breakpoints.


rxs.test <- data.frame(rx = c('dzyhapi', 'qqrvibe AB', 'aarvbyol', 'prlylyfe', 'ytrwow', 'zorncat',
					'qwexinox', 'xyxyslim', 'axyzoz XL', 'poyqwr'),
				loss = c(.04, .13, .02, .04, .06, .10, .09, .01, .05, .05))

exps.test <- data.frame(rx = c('dzyhapi', 'dzyhapi', 'dzyhapi', 'prlylyfe', 
					'qqrvibe AB', 'qqrvibe AB', 'qqrvibe AB', 
					'xyxyslim', 'xyxyslim', 'xyxyslim', 'xyxyslim', 
					'ytrwow', 'ytrwow'),
				Ns = c(15, 30, 36, 16, 26, 36, 42, 44, 20, 24, 19, 32, 29))



do.trickier.thing <- function(meds, exps, resample.n = 1000)
{
  # This does some sort of sampling-based metaanalytical computation. For each medicine,
  # there's several experiments, each with a different number of subjects. Do some
  # sort of silly computation that's a mean of those sample sizes plus the largest loss
  # from some other data. Lord knows why.
  
	# want those loss values
	merged <- merge(exps, meds, by='rx')

	# use resampling to get something probably meaningless
	meaning.less <- ddply(merged, .(rx), function(df) {
		data.frame(asdf = mean(sample(df$Ns, resample.n, replace=TRUE)) + max(df$loss))
	})

	meaning.less
}

do.trickier.thing(rxs.test, exps.test)



# 

piratical_dat.csv