This function performs the sample size calculation for a mixed model of
repeated measures with general correlation structure. See Lu, Luo, & Chen
(2008) for parameter definitions and other details. This function executes
Formula (3) on page 4.
Usage
power.mmrm(N = NULL, Ra = NULL, ra = NULL, sigmaa = NULL, Rb = NULL,
rb = NULL, sigmab = NULL, lambda = 1, delta = NULL,
sig.level = 0.05, power = NULL, alternative = c("two.sided",
"one.sided"), tol = .Machine$double.eps^2)
Arguments
N
total sample size
Ra
correlation matrix for group a
ra
retention in group a
sigmaa
standard deviation of observation of interest in group a
Rb
correlation matrix for group a
rb
retention in group b
sigmab
standard deviation of observation of interest in group b. If
NULL, sigmab is assumed same as sigmaa. If not NULL,
sigmaa and sigmab are averaged.
lambda
allocation ratio
delta
effect size
sig.level
type one error
power
power
alternative
one- or two-sided test
tol
numerical tolerance used in root finding.
Details
See Lu, Luo, & Chen (2008).
Value
The number of subject required per arm to attain the specified
power given sig.level and the other parameter estimates.
Author(s)
Michael C. Donohue
References
Lu, K., Luo, X., Chen, P.-Y. (2008) Sample size estimation for
repeated measures analysis in randomized clinical trials with missing data.
International Journal of Biostatistics, 4, (1)
See Also
power.mmrm.ar1, lmmpower,
diggle.linear.power
Examples
# reproduce Table 1 from Lu, Luo, & Chen (2008)
phi1 <- c(rep(1, 6), 2, 2)
phi2 <- c(1, 1, rep(2, 6))
lambda <- c(1, 2, sqrt(1/2), 1/2, 1, 2, 1, 2)
ztest <- ttest1 <- c()
for(i in 1:8){
Na <- (phi1[i] + lambda[i] * phi2[i])*(qnorm(0.05/2) + qnorm(1-0.90))^2*(0.5^-2)
Nb <- Na/lambda[i]
ztest <- c(ztest, Na + Nb)
v <- Na + Nb - 2
Na <- (phi1[i] + lambda[i] * phi2[i])*(qt(0.05/2, df = v) + qt(1-0.90, df = v))^2*(0.5^-2)
Nb <- Na/lambda[i]
ttest1 <- c(ttest1, Na + Nb)
}
data.frame(phi1, phi2, lambda, ztest, ttest1)
Ra <- matrix(0.25, nrow = 4, ncol = 4)
diag(Ra) <- 1
ra <- c(1, 0.90, 0.80, 0.70)
sigmaa <- 1
power.mmrm(Ra = Ra, ra = ra, sigmaa = sigmaa, delta = 0.5, power = 0.80)
power.mmrm(N = 174, Ra = Ra, ra = ra, sigmaa = sigmaa, delta = 0.5)
power.mmrm(N = 174, Ra = Ra, ra = ra, sigmaa = sigmaa, power = 0.80)
power.mmrm(Ra = Ra, ra = ra, sigmaa = sigmaa, delta = 0.5, power = 0.80, lambda = 2)
power.mmrm(N = 174, Ra = Ra, ra = ra, sigmaa = sigmaa, delta = 0.5, lambda = 2)
power.mmrm(N = 174, Ra = Ra, ra = ra, sigmaa = sigmaa, power = 0.80, lambda = 2)
# Extracting paramaters from gls objects with general correlation
# Create time index:
Orthodont$t.index <- as.numeric(factor(Orthodont$age, levels = c(8, 10, 12, 14)))
with(Orthodont, table(t.index, age))
fmOrth.corSym <- gls( distance ~ Sex * I(age - 11),
Orthodont,
correlation = corSymm(form = ~ t.index | Subject),
weights = varIdent(form = ~ 1 | age) )
summary(fmOrth.corSym)$tTable
C <- corMatrix(fmOrth.corSym$modelStruct$corStruct)[[1]]
sigmaa <- fmOrth.corSym$sigma *
coef(fmOrth.corSym$modelStruct$varStruct, unconstrained = FALSE)['14']
ra <- seq(1,0.80,length=nrow(C))
power.mmrm(N=100, Ra = C, ra = ra, sigmaa = sigmaa, power = 0.80)
# Extracting paramaters from gls objects with compound symmetric correlation
fmOrth.corCompSymm <- gls( distance ~ Sex * I(age - 11),
Orthodont,
correlation = corCompSymm(form = ~ t.index | Subject),
weights = varIdent(form = ~ 1 | age) )
summary(fmOrth.corCompSymm)$tTable
C <- corMatrix(fmOrth.corCompSymm$modelStruct$corStruct)[[1]]
sigmaa <- fmOrth.corCompSymm$sigma *
coef(fmOrth.corCompSymm$modelStruct$varStruct, unconstrained = FALSE)['14']
ra <- seq(1,0.80,length=nrow(C))
power.mmrm(N=100, Ra = C, ra = ra, sigmaa = sigmaa, power = 0.80)
# Extracting paramaters from gls objects with AR1 correlation
fmOrth.corAR1 <- gls( distance ~ Sex * I(age - 11),
Orthodont,
correlation = corAR1(form = ~ t.index | Subject),
weights = varIdent(form = ~ 1 | age) )
summary(fmOrth.corAR1)$tTable
C <- corMatrix(fmOrth.corAR1$modelStruct$corStruct)[[1]]
sigmaa <- fmOrth.corAR1$sigma *
coef(fmOrth.corAR1$modelStruct$varStruct, unconstrained = FALSE)['14']
ra <- seq(1,0.80,length=nrow(C))
power.mmrm(N=100, Ra = C, ra = ra, sigmaa = sigmaa, power = 0.80)
power.mmrm.ar1(N=100, rho = C[1,2], ra = ra, sigmaa = sigmaa, power = 0.80)