library(PracTools)
library(kableExtra)
# Model with quantitative covariates
data(PracTools::hospital)
HOSP <- PracTools::hospital
HOSP$sqrt.x <- sqrt(HOSP$x)
sam <- sample(nrow(HOSP), 50)
N1.resp <- HOSP[sam, ]
N2.nonresp <- HOSP[-sam, ]
## Create lm object using "known" data; no intercept model
lm.obj <- lm(y ~ 0 + sqrt.x + x, data = N1.resp)
del <- mean(HOSP$y) - mean(HOSP$y) * seq(.6, 1, by=0.05)
SampStop(lm.obj = lm.obj,
formula = ~ 0 + sqrt.x + x,
n1.data = N1.resp,
yvar = "y",
n2.data = N2.nonresp,
p = seq(0.2, 0.6, by=0.05),
delta = del,
seed = .Random.seed[413])
# Model with factors
data(PracTools::labor)
Labor <- PracTools::labor
sam <- sample(nrow(Labor), 50)
n1.vars <- c("WklyWage", "HoursPerWk", "agecat", "sex")
n2.vars <- c("HoursPerWk", "agecat", "sex")
N1.resp <- Labor[sam, n1.vars]
N2.nonresp <- Labor[-sam, n2.vars]
lm.obj <- lm(WklyWage ~ HoursPerWk + as.factor(agecat) + as.factor(sex), data = Labor)
del <- mean(N1.resp$WklyWage) - mean(N1.resp$WklyWage) * seq(.75, .95, by=0.05)
S <- SampStop(lm.obj = lm.obj,
formula = ~ HoursPerWk + as.factor(agecat) + as.factor(sex),
n1.data = N1.resp,
yvar = "WklyWage",
n2.data = N2.nonresp,
p = seq(0.2, 0.4, by=0.05),
delta = del,
seed = .Random.seed[78])
kableExtra::kable(S$Input, caption = "SampStop Input")
kableExtra::kable(head(S$Output, n=15),
caption = "SampStop Output: First 15 Observations")
S1 <- as.data.frame(S$Output)
p.nresp <- paste(S1$`Pr(response)`, S1$`Exp no. resps`, sep=", ")
library(ggplot2)
ggplot(S1, aes(x = `diff in means`,
y = `Pr(smaller diff)`,
colour = factor(p.nresp))) +
geom_point() +
geom_line(linewidth=1.1) +
labs(x = "delta", y = "Pr(|e1-e2|<= delta)", colour = "Pr(resp), n.resp")
Run the code above in your browser using DataLab