| Title: | Power and Sample Size for Linear Mixed Effect Models |
|---|---|
| Description: | Power and sample size calculation for testing fixed effect coefficients in multilevel linear mixed effect models with one or more than one independent populations. Laird, Nan M. and Ware, James H. (1982) <doi:10.2307/2529876>. |
| Authors: | Marco Chak Yan YU |
| Maintainer: | Marco Chak Yan YU <[email protected]> |
| License: | GPL-3 |
| Version: | 0.9.0 |
| Built: | 2026-05-12 06:11:08 UTC |
| Source: | https://github.com/cran/pass.lme |
Consider the following model:
Y=XB+Zu+e, u~N(0,D), e~N(0,R)
Yi~N(XBi,Vi), Vi=Zi*D*Zi'+Ri,
for each independent observation i
estimate of fixed effect coefficient B, denoted by b:
b=inv(sum(Xi'*inv(Vi)*Xi))*(sum(Xi'*inv(Vi)*Yi))
variance of b:
var(b)=Vb/n=inv(sum(Xi'*inv(Vi)*Xi))
where Vb=inv(Xi'*inv(Vi)*Xi)
distribution of any linear combinations L of b is given by:
Lb~N(mu,Sigma/n)
where mu = LB, Sigma = L*Vb*L'
lme.Lb.dist.theta(B, D, R, X, Z, m = NULL, L)lme.Lb.dist.theta(B, D, R, X, Z, m = NULL, L)
B |
fixed effect beta in px1 matrix |
D |
list of qxq random effect variance matrix;
where the first element corresponding to the highest-level effect,
the last element corresponding to the level 1 effect |
R |
residual variance |
X |
nxp matrix representing the covariates for the fixed effects |
Z |
nxq matrix representing the covariates for each level of random effects |
m |
vector of repeated measures from the highest to lowest level
(level 1 effects are addressed by Z and X and no need to be specified) |
L |
lxp matrix, representing l-linear-combinations of beta interested, |
theta: parameters (mu and Sigma) of the normal distribution for Lb
License: GPL-3
Marco Chak Yan YU
Maintainer: Marco Chak Yan YU <[email protected]>
#Example 1 # calc BLUE for 1-level LME model, # with covariates X, Z: (1,t), t=1,2,3 # for both fixed and random effects, # with fixed effect coefficients B: (100,-0.5), # random effect variance D: (2 1;1 2), # residual variance R: 0.2 B <- matrix(c(100,-0.5),2,1) D <- matrix(c(2,1,1,2),2,2) R <- 0.2 X <- cbind(rep(1,3),1:3) Z <- X lme.Lb.dist.theta(B,D,R,X,Z) #Example 2 # calc BLUE for 3-levels LME model, # with level 1 same as the above example # with 3 repeated-measures in level 2 # and 5 repeated-measures in the highest level, # assuming random effect variance for level 2 = (3 1;1 3), # and random effect variance for highest level = (5 1;1 5) D <- list(matrix(c(2,1,1,2),2,2),matrix(c(3,1,1,3),2,2), matrix(c(5,1,1,5),2,2)) lme.Lb.dist.theta(B,D,R,X,Z,m=c(5,3))#Example 1 # calc BLUE for 1-level LME model, # with covariates X, Z: (1,t), t=1,2,3 # for both fixed and random effects, # with fixed effect coefficients B: (100,-0.5), # random effect variance D: (2 1;1 2), # residual variance R: 0.2 B <- matrix(c(100,-0.5),2,1) D <- matrix(c(2,1,1,2),2,2) R <- 0.2 X <- cbind(rep(1,3),1:3) Z <- X lme.Lb.dist.theta(B,D,R,X,Z) #Example 2 # calc BLUE for 3-levels LME model, # with level 1 same as the above example # with 3 repeated-measures in level 2 # and 5 repeated-measures in the highest level, # assuming random effect variance for level 2 = (3 1;1 3), # and random effect variance for highest level = (5 1;1 5) D <- list(matrix(c(2,1,1,2),2,2),matrix(c(3,1,1,3),2,2), matrix(c(5,1,1,5),2,2)) lme.Lb.dist.theta(B,D,R,X,Z,m=c(5,3))
Interested parameters/linear combinations LB from more than one
independent populations can be aggregrate togeter by
appending mu vertically and Sigma/n diagonally
Consider Lb~N(MU,SIGMA) as the aggregrated estimates
Any comparison of interested parameters can be formulated by
multiplying a contrast matrix C on LB and set
H0: C*LB=d for any vector of value d to be tested
We then have
C*Lb~N(C*MU,C*SIGMA*C')
and
(C*Lb-d)'*inv(C*SIGMA*C')*(C*Lb-d)~chisq(q,lambda)
where degree of freedom q=rank(C*SIGMA*C'),
non-centrality parameter lambda=(C*LB-d)'*inv(C*SIGMA*C')*(C*LB-d)
Power of the test H0 is given by 1-beta=P(chisq(q,lambda)>qchisq(1-alpha,lambda))
Required sample size for desired power can be obtained by bisection method.
pass.lme.CLb.test(thetas, C = NULL, d = NULL, alpha = 0.05, power = NULL, n = NULL)pass.lme.CLb.test(thetas, C = NULL, d = NULL, alpha = 0.05, power = NULL, n = NULL)
thetas |
list of theta (LB and VLb), can be different for each group |
C |
Contrast of Matrix |
d |
Value vector to be tested for all contrast |
alpha |
significant level |
power |
desired power for sample size calculation |
n |
sample size for power calculation / |
solved.power given sample size n, this gives the power for testing H0
solved.n given the desired power, this gives the sample size for H0
License: GPL-3
Marco Chak Yan YU
Maintainer: Marco Chak Yan YU <[email protected]>
#Example 1 (test fixed effect coefficient 2=0) with power of 80% # for 1-level LME model, with covariates X, Z: (1,t), t=1,2,3 # for both fixed and random effects, with fixed effect coefficients B: (100,-0.5), # random effect variance D: (2 1;1 2), residual variance R: 0.2 B <- matrix(c(100,-0.5),2,1) D <- matrix(c(2,1,1,2),2,2) R <- 0.2 X <- cbind(rep(1,3),1:3) Z <- X theta <- lme.Lb.dist.theta(B,D,R,X,Z) pass.lme.CLb.test(list(theta),alpha=0.05,power=0.8) pass.lme.CLb.test(list(theta),alpha=0.05,n=66) #Example 2 (compare two fixed effect coefficient 2) with power of 80% # Consider above model as a control group model, # with an independent treatment group with model same as the control # except a different fixed effect coefficient 2 for treatment # = fixed effect coefficient 2 for control x 0.7 theta2 <- theta theta2$mu <- theta$mu *0.7 C <- matrix(c(1,-1),1,2) pass.lme.CLb.test(list(theta,theta2),C,alpha=0.05,power=0.8) pass.lme.CLb.test(list(theta,theta2),C,alpha=0.05,n=1468) #Example 3 (compare two fixed effect coefficient 2) with power of 80% # with sample size ratio, control:treatment = 1:2 pass.lme.CLb.test(list(theta,theta2),C,alpha=0.05,power=0.8,n=c(1,2)) pass.lme.CLb.test(list(theta,theta2),C,alpha=0.05,n=c(1101,2202)) #Example 4 (repeated-measures ANOVA for comparing 3 group means) with power of 80% # for 1-level LME model with mean for group 1, 2 and 3 are 100, 99, 102, respectively, # each subject to be measured 2 times, with within-subject variance = 15, residual variance = 10 B <- 100 D <- 15 R <- 10 X <- matrix(1,2,1) Z <- X theta <- lme.Lb.dist.theta(B,D,R,X,Z) theta2 <- theta theta3 <- theta theta2$mu <- 99 theta3$mu <- 102 C <- rbind(c(1,-1,0),c(1,0,-1)) pass.lme.CLb.test(list(theta,theta2,theta3),C,alpha=0.05,power=0.8) pass.lme.CLb.test(list(theta,theta2,theta3),C,alpha=0.05,n=41)#Example 1 (test fixed effect coefficient 2=0) with power of 80% # for 1-level LME model, with covariates X, Z: (1,t), t=1,2,3 # for both fixed and random effects, with fixed effect coefficients B: (100,-0.5), # random effect variance D: (2 1;1 2), residual variance R: 0.2 B <- matrix(c(100,-0.5),2,1) D <- matrix(c(2,1,1,2),2,2) R <- 0.2 X <- cbind(rep(1,3),1:3) Z <- X theta <- lme.Lb.dist.theta(B,D,R,X,Z) pass.lme.CLb.test(list(theta),alpha=0.05,power=0.8) pass.lme.CLb.test(list(theta),alpha=0.05,n=66) #Example 2 (compare two fixed effect coefficient 2) with power of 80% # Consider above model as a control group model, # with an independent treatment group with model same as the control # except a different fixed effect coefficient 2 for treatment # = fixed effect coefficient 2 for control x 0.7 theta2 <- theta theta2$mu <- theta$mu *0.7 C <- matrix(c(1,-1),1,2) pass.lme.CLb.test(list(theta,theta2),C,alpha=0.05,power=0.8) pass.lme.CLb.test(list(theta,theta2),C,alpha=0.05,n=1468) #Example 3 (compare two fixed effect coefficient 2) with power of 80% # with sample size ratio, control:treatment = 1:2 pass.lme.CLb.test(list(theta,theta2),C,alpha=0.05,power=0.8,n=c(1,2)) pass.lme.CLb.test(list(theta,theta2),C,alpha=0.05,n=c(1101,2202)) #Example 4 (repeated-measures ANOVA for comparing 3 group means) with power of 80% # for 1-level LME model with mean for group 1, 2 and 3 are 100, 99, 102, respectively, # each subject to be measured 2 times, with within-subject variance = 15, residual variance = 10 B <- 100 D <- 15 R <- 10 X <- matrix(1,2,1) Z <- X theta <- lme.Lb.dist.theta(B,D,R,X,Z) theta2 <- theta theta3 <- theta theta2$mu <- 99 theta3$mu <- 102 C <- rbind(c(1,-1,0),c(1,0,-1)) pass.lme.CLb.test(list(theta,theta2,theta3),C,alpha=0.05,power=0.8) pass.lme.CLb.test(list(theta,theta2,theta3),C,alpha=0.05,n=41)