t-test & Welch-test [Power / Sample Size]

Hi Rocco,

» Where does the formula on slide 10.83 in bebac.at/lectures/Leuven2013WS2.pdf for CI for parallel design come from? I cannot seem to find reference anywhere.

Honestly, I don’t remember why I simplified the commonly used formula.
Algebra:$$s\sqrt{\tfrac{n_1+n_2}{n_1n_2}}=\sqrt{s^2(1/n_1+1/n_2)}\;\tiny{\square}$$ Comparison with the data of the example.
• Like in the presentation:
mean.log <- function(x) mean(log(x), na.rm = TRUE)
T      <- c(100, 103, 80, 110,  78, 87, 116, 99, 122, 82, 68,  NA)
R      <- c(110, 113, 96,  90, 111, 68, 111, 93,  93, 82, 96, 137)
n1     <- sum(!is.na(T))
n2     <- sum(!is.na(R))
s1.2   <- var(log(T), na.rm = TRUE)
s2.2   <- var(log(R), na.rm = TRUE)
s0.2   <- ((n1 - 1) * s1.2 + (n2 - 1) * s2.2) / (n1 + n2 - 2)
s0     <- sqrt(s0.2)
nu.1   <- n1 + n2 - 2
nu.2   <- (s1.2 / n1 + s2.2 / n2)^2 /
(s1.2^2 / (n1^2 * (n1 - 1)) + s2.2^2 / (n2^2 * (n2 - 1)))
t.1    <- qt(p = 1 - 0.05, df = nu.1)
t.2    <- qt(p = 1 - 0.05, df = nu.2)
PE.log <- mean.log(T) - mean.log(R)
CI.1   <- PE.log + c(-1, +1) * t.1 *
s0 * sqrt((n1 + n2) / (n1 * n2))
CI.2   <- PE.log + c(-1, +1) * t.2 *
sqrt(s1.2 / n1 + s2.2 / n2)
CI.t   <- 100 * exp(CI.1)
CI.w   <- 100 * exp(CI.2)
fmt    <- "%.3f %.3f %.3f %.2f    %.2f   %.2f"
cat("    method     df mean.T mean.R    PE CL.lower CL.upper",
"\n    t-test",
sprintf(fmt, nu.1, exp(mean.log(T)), exp(mean.log(R)),
100 * exp(PE.log), CI.t, CI.t),
"\nWelch-test",
sprintf(fmt, nu.2, exp(mean.log(T)), exp(mean.log(R)),
100 * exp(PE.log), CI.w, CI.w), "\n")

method     df mean.T mean.R    PE CL.lower CL.upper
t-test 21.000 93.554 98.551 94.93    83.28   108.20
Welch-test 20.705 93.554 98.551 94.93    83.26   108.23

• More comfortable with the t.test() function, where var.equal = TRUE gives the t-test and var.equal = FALSE the Welch-test:
res <- data.frame(method = c("t-test", "Welch-test"), df = NA,
mean.T = NA, mean.R = NA, PE = NA,
CL.lower = NA, CL.upper = NA)
var.equal <- c(TRUE, FALSE)
for (j in 1:2) {
x <- t.test(x = log(T), y = log(R), conf.level = 0.90,
var.equal = var.equal[j])
res[j, 2]   <- signif(x[], 5)
res[j, 3:4] <- signif(exp(x[][1:2]), 5)
res[j, 5]   <- round(100 * exp(diff(x[][2:1])), 2)
res[j, 6:7] <- round(100 * exp(x[]), 2)
}
print(res, row.names = FALSE)

method     df mean.T mean.R    PE CL.lower CL.upper
t-test 21.000 93.554 98.551 94.93    83.28   108.20
Welch-test 20.705 93.554 98.551 94.93    83.26   108.23
The formula for Satterthwaite’s approximation of the degrees of freedom given in slide 11 contained typos (corrected in the meantime). Of course,$$\nu\approx\frac{\left(\frac{{s_{1}}^{2}}{n_1}+\frac{{s_{2}}^{2}}{n_2}\right)^2}{\frac{{s_{1}}^{4}}{n{_{1}}^{2}(n_1-1)}+\frac{{s_{2}}^{4}}{n{_{2}}^{2}(n_2-1)}}$$Satterthwaite’s approximation adjusts both for unequal variances and group sizes. The conventional t-test is fairly robust against the former but less so for the latter.
In the R function t.test() var.equal = FALSE is the default because:
• Any pre-test is bad practice (might inflate the type I error) and has low power esp. for small sample sizes.
• If ${s_{1}}^{2}={s_{2}}^{2}\;\wedge\;n_1=n_2$, the formula given above reduces to the simple $\nu=n_1+n_2-2$ anyhow.
• In all other cases the Welch-test is conservative, which is a desirable property.

Cheers,
Helmut Schütz The quality of responses received is directly proportional to the quality of the question asked. ☼
Science Quotes

Bioequivalence and Bioavailability Forum |  Admin contact
19,871 posts in 4,212 threads, 1,364 registered users;
online 8 (0 registered, 8 guests [including 7 identified bots]).
Forum time (Europe/Vienna): 12:40 CEST

A drug is not bad. A drug is a chemical compound.
The problem comes in when people who take drugs treat them
like a license to behave like an asshole.    Frank Zappa

The BIOEQUIVALENCE / BIOAVAILABILITY FORUM is hosted by Ing. Helmut Schütz 