Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Code example for post `Randomization-based inference: the Neymanian approach'
library(magrittr)
library(purrr)
library(gtools)
## Estimand function and estimator function
tau.fn <- function(science) mean(science$y1) - mean(science$y0)
hat.tau.fn <- function(w, y) mean(y[w==1]) - mean(y[w==0])
## ----------------------------
## The science table + estimand
## ----------------------------
science <- data.frame(y0 = c(2, 1, 2, 0, 3, 0),
y1 = c(3, 1, 4, 1, 3, 2))
tau <- tau.fn(science)
## ----------------------------
## observed data
## ----------------------------
W.obs <- c(1, 1, 0, 1, 0, 0)
Y.obs <- science$y1 * W.obs + science$y0 * (1-W.obs)
hat.tau.obs <- hat.tau.fn(W.obs, Y.obs)
## ----------------------------
## bias of estimator
## ----------------------------
W.ls <- permutations(6, 6, 1:6) %>% array_tree(1) %>% map(~ W.obs[.]) %>% unique
hat.tau.ls <- vector('numeric', length=length(W.ls))
for(i in seq_along(W.ls)) {
W <- W.ls[[i]]
Y <- science$y1 * W + science$y0 * (1-W)
hat.tau.ls[i] <- hat.tau.fn(W, Y)
print(hat.tau.fn(W, Y))
}
bias <- mean(hat.tau.ls) - tau
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment