swGlmPwr | R Documentation |
swGlmPwr
returns (two-sided) power of the treatment effect for the specified SW CRT design in the context of generalized linear models by adopting the Laplace approximation detailed in Breslow and Clayton (1993) to obtain the covariance matrix of the estimated parameters. The response/outcome of interest can be binomial (logit link only) or Poisson (log link only) distributed. A cross-sectional or closed cohort sampling scheme can be specified.
The outcome is assumed to come from a model with fixed treatment effect(s) (using an immediate treatment (IT) or exposure time indicator (ETI) model - see Kenny et al (2022)), fixed time effect, random intercepts, random treatment effects, random cluster-specific time effects and, in the case of closed cohort sampling, an individual random effect. The coefficients for fixed effects can be specified using fixed.intercept
, fixed.treatment.effect
, and fixed.time.effect
. Variance components can be specified using tau
, eta
, rho
, gamma
and zeta
.
swGlmPwr(design, distn, n, fixed.intercept, fixed.treatment.effect,
fixed.time.effect, H = NULL, tau = 0, eta = 0, rho = 0, gamma = 0, zeta = 0, ar=1,
alpha=0.05, retDATA = FALSE, silent=FALSE)
design |
list: A stepped wedge design object, typically from |
distn |
character: Distribution assumed (binomial or Poisson). "binomial" implies binomial outcomes with logit link and "poisson" implies Poisson outcome with log link. ***NOTE: It is the users responsibility to make sure specified parameters (fixed.intercept, fixed.treatment effect, fixed.time effect, tau, eta, rho, gamma, zeta) are ALL on SAME scale as specified link function; see example.*** |
n |
integer (scalar, vector, or matrix): Number of observations: (scalar) for all clusters and all time points that are not NA in design$swDsn; (vector) for each cluster at all time points that are not NA in design$swDsn; and (matrix) for each cluster at each time point, where rows correspond to clusters and columns correspond to time. |
fixed.intercept |
numeric (scalar): Intercept for the fixed effect on canonical scales (logit for binomial outcomes and log for Poisson outcomes). It is the mean outcome under the control condition in the first time point transformed to the canonical scales. |
fixed.treatment.effect |
numeric (scalar, vector): Gives the coefficients for the treatment effect on the canonical scale (logit for binomial outcomes and log for Poisson outcomes). If If |
fixed.time.effect |
numeric(scalar, vector): Coefficients for the time (as dummy variables) in the fixed effect model on canonical scales (logit for binomial outcomes and log for Poisson outcomes). The first time point is always used as reference. Specify a common time effect for all time points after the first (scalar) or differnt time effects for each time point after the first (vector of length (total time-1)). |
H |
numeric (vector): If NULL, then swGlmPwr assumes an immediate, constant treatment effect (IT) model. If not NULL, then an exposure time indicator (ETI) model is assumed and H should be a vector as long as the longest exposure time (in a classic SW design, the longest exposure period is the number of time periods minus one). H specifies the desired linear combination of fixed.treatment.effect. For example, in a stepped wedge trial with 5 time periods and four exposure times, H = rep(.25,4) gives the average treatment effect over the four exposure times; H = c(0,0,.5,.5) ignores the first two periods after the intervention is introduced and averages the remaining periods. Typically, the sum of H is 1.0; if not, it is renormalized to sum to 1.0. H can only be specified when there is a single intervention type (i.e. design$swDsn includes only NA,0,1; see cautions in help for |
tau |
numeric (scalar): Standard deviation of random intercepts on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
eta |
numeric (scalar): Standard deviation of random treatment effects on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
rho |
numeric (scalar): Correlation between random intercepts and random treatment effects on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
gamma |
numeric (scalar): Standard deviation of random time effects on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
zeta |
numeric (scalar): Standard deviation of individual effects for closed cohort sampling on canonical scales (logit for binomial outcomes and log for Poisson outcomes). Default is 0 which implies cross-sectional sampling. Values greater than 0 imply closed cohort sampling. |
ar |
numeric (scalar): Autocorrelation for exponential decay structure between periods within a cluster as described in Kasza et al (2019a;2019b). Available for cross-sectional sampling only (i.e. zeta=0 and iac=0). Default is ar=1 (no exponential decay). See details. |
alpha |
numeric (scalar): Statistical significance level. Default is 0.05. |
retDATA |
logical: if |
silent |
logical: if TRUE, hides most warning messages. Default value is |
Let Y_i
denote a vector of observations for cluster i
. The generalized linear mixed model formulation assumes
E(Y_i | b_i) = h(X\beta + Z_ib_i)
where h
is a link function, \beta
is a vector of fixed effects, b_i
is a vector of random effects (in order, cluster, cluster*time, cluster*treatment) with b_i \sim N(0,D)
and X
and Z
are design matricies for the fixed and random effects, respectively. Then Breslow and Clayton (1993) show that
Var(\hat{\beta}) \approx \left( X^T V^{-1} X \right)^{-1}
where
V = W^b + ZDZ^T
where W^b
is a diagonal matrix (W^b
is approximated as W^0
in all calculations).
An exponentially decaying autocorrelation structure (only available for cross-sectional sampling i.e. zeta
= 0) is incorporated by modifying V
to
V = (W^b + ZD1Z^T) * R + ZD2Z^T
where *
denotes element-by-element matrix multiplication and, assuming there are random effects for cluster, cluster*time and cluster*treatment,
D1_{ij} = D
for i and j <= J+1 (corresponding to the cluster and cluster*time random effects) and 0 otherwise
D2_{ij} = D
for i or j > J+1 (corresponding to the cluster*treatment random effect) and 0 otherwise
R
= a toeplitz matrix with elements equal to ar
^m where m
= 0...J-1 is the distance from the main diagonal.
The two-sided statistical power for treatment effect \theta
(equal to H
%*%fixed.treatment.effect
if H
is non-NULL) is
Pwr(\theta) = \Phi( \frac{Z - z_{1 - \alpha /2} \sqrt{V_0(\hat{\theta})}}{\sqrt{V_\alpha(\hat{\theta})}}) + 1 - \Phi( \frac{Z+ z_{1 - \alpha /2} \sqrt{V_0(\hat{\theta})}}{\sqrt{V_\alpha(\hat{\theta})}})
,
where \Phi
is the cumulative distribution function of the standard normal distribution.
The variance of \hat{\theta}
under the null is denoted as V_0(\hat{\theta})
, and the variance of \hat{\theta}
under the alternative is denoted as V_\alpha(\hat{\theta})
). Both variances are approximated by simplifying the Laplace approximation that marginalizes the random effects in the generalized linear mixed models. For more details, see Xia et al. (2020).
When the outcome is Gaussian, the method adopted by swGlmPwr
coincides with that of swPwr
, so power calculation for Gaussian outcomes is not included in swGlmPwr
to avoid repetition. When the outcome is binomial, swGlmPwr
performs power calculation on the natural scale (logit), while swPwr
performs power calculation on the linear scale.
The value of zeta defines the samling scheme. When zeta = 0, cross-sectional sampling is assumed; if zeta > 0 then closed cohort sampling is assumed.
Incomplete designs (e.g. staircase) can be specified by using the swDsn
function with the swBlk
argument to specify the incomplete design or by specifying a complete design with swDsn
and using the n
argument of swPwr to specify a matrix with 0 in the missing cluster-periods. Note that if the second approach is used then H must be either a scalar or a vector as long as the longest lag in the complete design. See the examples below.
numeric (scalar): swGlmPwr
returns the power of treatment effect if retDATA = FALSE.
numeric (list): swGlmPwr
returns all specified and computed items as objects of a list if retDATA = TRUE.
$design |
list: A stepped wedge design object, typically from swDsn, that includes at least the following components: swDsn, swDsn.unique.clusters, clusters, n.clusters, total.time, nTxLev |
$distn |
character: Distribution assumed (binomial or Poisson). "binomial" implies binomial outcomes and "poisson" implies Poisson outcome. |
$n |
integer (scalar, vector, or matrix): Number of observations: (scalar) for all clusters and all time points; (vector) for each cluster at all time points; and (matrix) for each cluster at each time point, where rows correspond to clusters and columns correspond to time. |
$fixed.intercept |
numeric (scalar): Intercept for the fixed effect on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
$fixed.treatment.effect |
numeric (scalar): Coefficient(s) for the treatment(s) in the fixed effect model on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
$fixed.time.effect |
numeric(scalar, vector): Coefficients for the time (as dummy variables) in the fixed effect model on canonical scales (logit for binomial outcomes and log for Poisson outcomes). The first time point is always used as reference. A common time effect for all time points after the first (scalar) or differnt time effects for each time point after the first (vector of length (total time-1)). |
H |
numeric (vector): H specifies the desired linear combination of exposure time treatment effects for a ETI model-based estimate. |
$tau |
numeric (scalar): Standard deviation of random intercepts on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
$eta |
numeric (scalar): Standard deviation of random treatment effects on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
$rho |
numeric (scalar): Correlation between random intercepts and random treatment effects on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
$gamma |
numeric (scalar): Standard deviation of random time effects on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
$zeta |
numeric (scalar): Standard deviation of individual random effects for closed cohort sampling on canonical scales (logit for binomial outcomes and log for Poisson outcomes). |
ar |
numeric (scalar): Autocorrelation for exponential decay structure between periods within a cluster as described in Kasza et al (2019). Available for cross-sectional sampling only (i.e. zeta=0 and iac=0). See details. |
$alpha |
numeric (scalar): Statistical significance level. Default is 0.05. |
$var.theta.null |
numeric (martix): Variance-covariance martix of the estimated treatment effect(s) under the null for this SW CRT design. |
$var.theta.alt |
numeric (marix): Variance-covariance matrix of the estimated treatment effect(s) under the alternative for this SW CRT design. |
$pwrGLM |
numeric (scalar): Power of treatment effect using a simplified Laplace approximation. |
Fan Xia, James P Hughes, and Emily C Voldal
Breslow NE and Clayton DG (1993). Approximate inference in generalized linear mixed models. Journal of the American Statistical Association, 88(421):9-25.
Kasza J, Hemming K, Hooper R, Matthews JNS, Forbes AB. Impact of non-uniform correlation on sample size and power in multiple-period cluster randomized trials. Statistical Methods in Medical Research 2019a; 28: 703-716.
Kasza J, Taljaard M, Forbes AB. Information content of stepped-wedge designs when treatment effect heterogeneity and/or implementation periods are present. Statistics in Medicine 2019b; 38:4686-4701.
Kenny A, Voldal E, Xia F, Heagerty PJ, Hughes JP. Analysis of stepped wedge cluster randomized trials in the presence of a time-varying treatment effect. Statistics in Medicine, in press, 2022.
Voldal EC, Hakhu NR, Xia, F, Heagerty PJ, Hughes JP. swCRTdesign: An R Package for Stepped Wedge Trial Design and Analysis. Computer Methods and Programs in Biomedicine 2020;196:105514.
Xia F, Hughes JP, Voldal EC, Heagerty PJ. Power and sample size calculation for stepped-wedge designs with discrete outcomes. Trials. 2021 Dec;22(1):598.
##test-case large clusters
logit <- function(x){log(x/(1 - x))}
#specify large cluster sizes
size = c(35219,53535,63785,456132,128670,96673,
51454,156667,127440,68615,56502,17719,75931,58655,52874,75936)
# Cross-sectional example
swGlmPwr(design=swDsn(c(4,3,5,4)),distn="binomial",n=size,
fixed.intercept=log(28.62/(2*100000)),fixed.time.effect = 1,
fixed.treatment.effect = log(.6),
tau=.31,eta=abs(0.4*log(.6)),rho=0,gamma=.15,alpha=.05,retDATA = FALSE)
# Closed cohort example, comparing average of intervention period exposure
# lags 3 and 4 to control period
swGlmPwr(design=swDsn(c(5,5,5,5,5),extra.ctl.time=3,extra.trt.time=5),
distn="binomial",n=20,
fixed.intercept=logit(0.40),
fixed.treatment.effect=c(0,0,rep(logit(0.60)-logit(0.40),8)),
fixed.time.effect=0.08,
H=c(0,0,1,1,0,0,0,0,0,0),
tau=sqrt(.1316),gamma=sqrt(.1974),eta=0,zeta=sqrt(2.5))
# Example wih periods with no data and multiple treatment levels
stdy <- swDsn(c(6,6,6,6),swBlk=matrix(c(0,1,2,2,2,2,
NA,0,1,2,2,2,
NA,NA,0,1,2,2,
NA,NA,NA,0,1,2),4,6,byrow=TRUE))
swGlmPwr(stdy, distn="binomial",n=10000,
fixed.intercept=log(28.62/(2*100000)),fixed.time.effect = 1,
fixed.treatment.effect = c(log(.6),log(.5)),
tau=.31,eta=0,rho=0,gamma=0,alpha=.05,retDATA = FALSE)
#Incomplete design
# These all give the same power
design12 = swDsn(c(4,4,4),swBlk=matrix(c(0,1,1,NA,NA,NA,0,1,1,NA,NA,NA,0,1,1),3,5,byrow=TRUE))
swGlmPwr(design12,distn="binomial",n=10,fixed.intercept=0,fixed.treatment.effect=c(.8,.8),
fixed.time.effect=0,H=c(.5,.5),tau=.1)
swGlmPwr(design12,distn="binomial",n=10,fixed.intercept=0,fixed.treatment.effect=c(.8,.8),
fixed.time.effect=0,H=1,tau=.1) # H is automatically expanded and normalized
#
design12a = swDsn(c(4,4,4),extra.trt.time=1)
nmat12a=matrix(c(rep(c(10,10,10,0,0),4),rep(c(0,10,10,10,0),4),rep(c(0,0,10,10,10),4)),12,5,
byrow=TRUE)
#Note maximum exposure time is 4 by design even though there are no data for exposure times 3 and 4
swGlmPwr(design12a,distn="binomial",n=nmat12a,fixed.intercept=0,
fixed.treatment.effect=c(.8,.8,.8,.8),fixed.time.effect=0,H=1,tau=.1) #OK
swGlmPwr(design12a,distn="binomial",n=nmat12a,fixed.intercept=0,
fixed.treatment.effect=c(.8,.8,.8,.8),fixed.time.effect=0,H=c(1,1,0,0),tau=.1) #OK
# The following generates an error; since design12a is a complete design with maximum lag 4,
# fixed.treatment.effect must be a vector length 4 and H must be a scalar or vector of length 4
# swGlmPwr(design12a,distn="binomial",n=nmat12a,fixed.intercept=0,fixed.treatment.effect=c(.8,.8),
# fixed.time.effect=0,H=c(1,1),tau=.1)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.