seed_value <- 19710822 set.seed(seed_value) barplot(dbinom(0:20, size = 20, prob = 0.5), names.arg = 0:20) # The null hypothesis is "prob = 0.5" 1 - (pbinom(12, 20, p = 0.5 ) - pbinom(7, 20, p = 0.5)) binom.test(13, 20, p = 0.5) binom.test(130, 200) x <- rbinom(10000, 20, p = 0.5) pval <- vector(length = 10000) for (n in 1:10000) { pval[n] <- binom.test(x[n], 20, p = 0.5)$p.value } print(sum(pval <= 0.05) / length(x)) # sum of booleans is number of TRUEs x <- seq(-4, 4, by = 0.01) plot(x, dnorm(x), type = "l", ylab="density") lines(x, dt(x, df = 4), type = "l", col = "red", lty=2) lines(x, dt(x, df = 10), type = "l", col = "blue", lty=4) legend(x = 1.25, y = 0.4, legend = c("normal", "df = 4", "df = 10"), lty = c(1, 4, 5), col = c("black", "red", "blue")) x1 <- rnorm(100); x2 <- rnorm(100); t.test(x2, x1) x1 <- rnorm(100); x2 <- rnorm(100, mean = 1); t.test(x2, x1) x1 <- rnorm(100); x2 <- rnorm(100, mean = 0.2); t.test(x2, x1) set.seed(2) n <- 3 x <- rnorm(n) y <- rnorm(n) pval <- c(NA, NA, t.test(x, y)$p.value) while(pval[n] > 0.05) { x[n + 1] <- rnorm(1) y[n + 1] <- rnorm(1) test <- t.test(x, y) pval[n + 1] <- test$p.value n <- n + 1 } plot(1:n, pval)