Playing with plotting difference-of-means tests

explot.R

set.seed(1)
n <- 300
x <- rbinom(n, 1, .5)
y <- 4 + (2 * x) + rnorm(n, 0, 3)
m <- tapply(y, x, mean)

wm <- weighted.mean(c(mean(y[x==1]), mean(y[x==0])), c(length(y[x==1]), length(y[x==0])))
ci <- function(formula, conf.level, var.equal = FALSE) {
    (tt <- t.test(formula, var.equal = var.equal, conf.level = conf.level))
    rect(0, wm - diff(tt$conf.int)/2,
         1, wm + diff(tt$conf.int)/2, col = rgb(1,0,0,.05), border = NA)
}
plot(x,y, col = gray(.5, .2), xlim=c(-1,2), ylim = c(2,8), xaxt = "n", las = 1)
axis(1, 0:1, 0:1)
c0 <- (sd(y[x==0])/sqrt(length(x[x==0])))
c1 <- (sd(y[x==1])/sqrt(length(x[x==1])))
segments(0, m[1] - 1 * c0, 0, m[1] + 1 * c0, lwd = 3, col = "blue")
segments(1, m[2] - 1 * c1, 1, m[2] + 1 * c1, lwd = 3, col = "blue")
segments(0, m[1] - 2 * c0, 0, m[1] + 2 * c0, lwd = 1, col = "blue")
segments(1, m[2] - 2 * c1, 1, m[2] + 2 * c1, lwd = 1, col = "blue")
points(0:1, m, col = "red", pch = 19, cex = 1)
invisible(lapply(seq(.1,.99, length.out = 20), ci, formula = y~x))
segments(0, wm, 1, wm, lwd = 1.5)