Description Usage Arguments Details Value Author(s) References Examples
Power calculation for simple linear mixed effects model. This function can be used to calculate one of the 3 parameters (power, sample size, and minimum detectable slope) by setting the corresponding parameter as NULL and providing values for the other 2 parameters.
1 2 3 4 5 6 7 8 9 10 11 12 |
slope |
numeric. Slope under alternative hypothesis. |
n |
integer. Total number of subjects. |
m |
integer. Number of observations per subject. |
sigma.y |
numeric. Standard deviation of the outcome y. |
sigma.x |
numeric. Standard deviation of the predictor x. |
power |
numeric. Desired power. |
rho |
numeric. Intra-class correlation (i.e., correlation between y_{ij} and y_{ik} for the j-th and k-th observations of the i-th subject). |
FWER |
numeric. Family-wise Type I error rate. |
nTests |
integer. Number of tests (e.g., number of genes in differential expression analysis based on scRNAseq to compare gene expression between diseased subjects and healthy subjects). |
n.lower |
numeric. Lower bound of the total number of subjects. Only used when calculating toal number of subjects. |
n.upper |
numeric. Upper bound of the total number of subjects. Only used when calculating total number of subjects. |
We assume the following simple linear mixed effects model to characterize the association between the predictor x and the outcome y:
y_{ij} = β_{0i} + β_1 * x_i + ε_{ij},
where
β_{0i} \sim N≤ft(β_0, σ^2_{β}\right),
and
ε_{ij} \sim N≤ft(0, σ^2\right),
i=1,…, n, j=1,…, m, n is the number of subjects, m is the number of observations per subject, y_{ij} is the outcome value for the j-th observation of the i-th subject, x_i is the predictor value for the i-th subject. For example, x_i is the binary variable indicating if the i-th subject is a diseased subject or not.
We would like to test the following hypotheses:
H_0: β_1=0,
and
H_1: β_1 = δ,
where δ\neq 0.
We can derive the power calculation formula is
power=1- Φ≤ft(z_{α^{*}/2}-a\times b\right) +Φ≤ft(-z_{α^{*}/2} - a\times b\right),
where
a= \frac{\hat{σ}_x }{σ_y}
and
b=\frac{δ√{m(n-1)}}{√{1+(m-1)ρ}}
and z_{α^{*}/2} is the upper 100α^{*}/2 percentile of the standard normal distribution, α^{*}=α/nTests, nTests is the number of tests, σ_y=√{σ^2_{β}+σ^2}, \hat{σ}_x=√{∑_{i=1}^n≤ft(x_i-\bar{x}\right)^2/(n-1)}, and ρ=σ^2_{β}/≤ft(σ^2_{β}+σ^2\right) is the intra-class correlation.
power if the input parameter power = NULL
.
sample size (total number of subjects) if the input parameter n = NULL
;
minimum detectable slope if the input parameter slope = NULL
.
Xianjun Dong <XDONG@rics.bwh.harvard.edu>, Xiaoqi Li<xli85@bwh.harvard.edu>, Tzuu-Wang Chang <Chang.Tzuu-Wang@mgh.harvard.edu>, Scott T. Weiss <restw@channing.harvard.edu>, Weiliang Qiu <weiliang.qiu@gmail.com>
Dong X, Li X, Chang T-W, Scherzer CR, Weiss ST, and Qiu W. powerEQTL: An R package and shiny application for sample size and power calculation of bulk tissue and single-cell eQTL analysis. Bioinformatics, 2021;, btab385
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 | n = 102
m = 227868
# calculate power
power = powerLME(
slope = 0.6,
n = n,
m = m,
sigma.y = 0.29,
sigma.x = 0.308,
power = NULL,
rho = 0.8,
FWER = 0.05,
nTests = 1e+6)
print(power)
# calculate sample size (total number of subjects)
n = powerLME(
slope = 0.6,
n = NULL,
m = m,
sigma.y = 0.29,
sigma.x = 0.308,
power = 0.9562555,
rho = 0.8,
FWER = 0.05,
nTests = 1e+6)
print(n)
# calculate slope
slope = powerLME(
slope = NULL,
n = n,
m = m,
sigma.y = 0.29,
sigma.x = 0.308,
power = 0.9562555,
rho = 0.8,
FWER = 0.05,
nTests = 1e+6)
print(slope)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.