liu.liang.linear.power: Linear mixed model sample size calculations from Liu & Liang (1997).

Description

This function performs the sample size calculation for a linear mixed model. See Liu and Liang (1997) for parameter definitions and other details.

Usage

liu.liang.linear.power(
  N = NULL,
  delta = NULL,
  u = NULL,
  v = NULL,
  sigma2 = 1,
  R = NULL,
  R.list = NULL,
  sig.level = 0.05,
  power = NULL,
  Pi = rep(1/length(u), length(u)),
  alternative = c("two.sided", "one.sided"),
  tol = .Machine$double.eps^2
)

Arguments

N: The total sample size. This formula can accommodate unbalanced group allocation via Pi. See Liu and Liang (1997) for more details
delta: group difference (possibly a vector of differences)
u: a list of covariate vectors or matrices associated with the parameter of interest
v: a respective list of covariate vectors or matrices associated with the nuisance parameter
sigma2: the error variance
R: the variance-covariance matrix for the repeated measures
R.list: a list of variance-covariance matrices for the repeated measures, if assumed different in two groups
sig.level: type one error
power: power
Pi: the proportion of covariates of each type
alternative: one- or two-sided test
tol: numerical tolerance used in root finding.

Details

The parameters u, v, and Pi are expected to be the same length and sorted with respect to each other. See Liu and Liang (1997) and package vignette for more details.

References

Liu, G. and Liang, K. Y. (1997) Sample size calculations for studies with correlated observations. Biometrics, 53(3), 937-47.

Examples

Run this code


if (FALSE) {
browseVignettes(package = "longpower")
}

# Reproduces the table on page 29 of Diggle et al for
# difference in slopes between groups

n <- 3
t <- c(0,2,5)
u <- list(u1 = t, u2 = rep(0,n))
v <- list(v1 = cbind(1,1,t),
         v2 = cbind(1,0,t))         
rho <- c(0.2, 0.5, 0.8)
sigma2 <- c(100, 200, 300)
tab <- outer(rho, sigma2, 
      Vectorize(function(rho, sigma2){
        ceiling(liu.liang.linear.power(
          delta=0.5, u=u, v=v,
          sigma2=sigma2,
          R=rho, alternative="one.sided",
          power=0.80)$N/2)}))
colnames(tab) <- paste("sigma2 =", sigma2)
rownames(tab) <- paste("rho =", rho)
tab

# Reproduces the table on page 30 of Diggle et al for 
# difference in average response between groups.

n <- 3
u <- list(u1 = rep(1,n), u2 = rep(0,n))
v <- list(v1 = rep(1,n),
         v2 = rep(1,n))
rho <- c(0.2, 0.5, 0.8)
delta <- c(20, 30, 40, 50)/100
tab <- outer(rho, delta, 
     Vectorize(function(rho, delta){
       ceiling(liu.liang.linear.power(
         delta=delta, u=u, v=v,
         sigma2=1,
         R=rho, alternative="one.sided",
         power=0.80)$n[1])}))
colnames(tab) <- paste("delta =", delta)
rownames(tab) <- paste("rho =", rho)
tab

# An Alzheimer's Disease example using ADAS-cog pilot estimates
# var of random intercept
sig2.i <- 55
# var of random slope
sig2.s <- 24
# residual var
sig2.e <- 10
# covariance of slope and intercep
cov.s.i <- 0.8*sqrt(sig2.i)*sqrt(sig2.s)

cov.t <- function(t1, t2, sig2.i, sig2.s, cov.s.i){
        sig2.i + t1*t2*sig2.s + (t1+t2)*cov.s.i 
}

t <- seq(0,1.5,0.25)
n <- length(t)
R <- outer(t, t, function(x,y){cov.t(x,y, sig2.i, sig2.s, cov.s.i)})
R <- R + diag(sig2.e, n, n)
u <- list(u1 = t, u2 = rep(0,n))
v <- list(v1 = cbind(1,1,t),
         v2 = cbind(1,0,t))         

liu.liang.linear.power(delta=1.5, u=u, v=v, R=R, sig.level=0.05, power=0.80)
liu.liang.linear.power(N=416, u=u, v=v, R=R, sig.level=0.05, power=0.80)
liu.liang.linear.power(N=416, delta = 1.5, u=u, v=v, R=R, sig.level=0.05)
liu.liang.linear.power(N=416, delta = 1.5, u=u, v=v, R=R, power=0.80, sig.level = NULL)

# Reproduces total sample sizes, m, of Table 1 of Liu and Liang 1997
tab1 <- data.frame(cbind(
  n = c(rep(4, 4), rep(2, 4), 1),
  rho = c(0.0, 0.3, 0.5, 0.8)))
m <- c()
for(i in 1:nrow(tab1)){
  R <- matrix(tab1$rho[i], nrow = tab1$n[i], ncol = tab1$n[i])
  diag(R) <- 1
  m <- c(m, ceiling(liu.liang.linear.power(
    delta=0.5,
    u = list(u1 = rep(1, tab1$n[i]), # treatment
             u2 = rep(0, tab1$n[i])), # control       
    v = list(v1 = rep(1, tab1$n[i]), v2 = rep(1, tab1$n[i])), # intercept
    sigma2=1,
    R=R, alternative="two.sided",
    power=0.90)$N))
}
cbind(tab1, m)

# Reproduces total sample sizes, m, of Table 3.a. of Liu and Liang 1997
# with unbalanced design
tab3 <- data.frame(cbind(
  rho = rep(c(0.0, 0.3, 0.5, 0.8), 2),
  pi1 = c(rep(0.8, 4), rep(0.2, 4))))
m <- c()
for(i in 1:nrow(tab3)){
  R <- matrix(tab3$rho[i], nrow = 4, ncol = 4)
  diag(R) <- 1
  m <- c(m, ceiling(liu.liang.linear.power(
    delta=0.5,
    u = list(u1 = rep(1, 4), # treatment
             u2 = rep(0, 4)), # control       
    v = list(v1 = rep(1, 4), v2 = rep(1, 4)), # intercept
    sigma2=1,
    Pi = c(tab3$pi1[i], 1-tab3$pi1[i]),
    R=R, alternative="two.sided",
    power=0.90)$N))
}
cbind(tab3, m)