powerLME {powerEQTL} | R Documentation |
Power Calculation for Simple Linear Mixed Effects Model
Description
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.
Usage
powerLME(
slope,
n,
m,
sigma.y,
sigma.x,
power = NULL,
rho = 0.8,
FWER = 0.05,
nTests = 1,
n.lower = 2.01,
n.upper = 1e+30)
Arguments
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
|
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. |
Details
We assume the following simple linear mixed effects model to characterize the association between the predictor x and the outcome y:
y_{ij} = \beta_{0i} + \beta_1 * x_i + \epsilon_{ij},
where
\beta_{0i} \sim N\left(\beta_0, \sigma^2_{\beta}\right),
and
\epsilon_{ij} \sim N\left(0, \sigma^2\right),
i=1,\ldots, n
, j=1,\ldots, 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: \beta_1=0,
and
H_1: \beta_1 = \delta,
where \delta\neq 0
.
We can derive the power calculation formula is
power=1- \Phi\left(z_{\alpha^{*}/2}-a\times b\right)
+\Phi\left(-z_{\alpha^{*}/2} - a\times b\right),
where
a=
\frac{\hat{\sigma}_x
}{\sigma_y}
and
b=\frac{\delta\sqrt{m(n-1)}}{\sqrt{1+(m-1)\rho}}
and z_{\alpha^{*}/2}
is the upper 100\alpha^{*}/2
percentile of the standard normal distribution,
\alpha^{*}=\alpha/nTests
, nTests is the number of
tests,
\sigma_y=\sqrt{\sigma^2_{\beta}+\sigma^2}
,
\hat{\sigma}_x=\sqrt{\sum_{i=1}^n\left(x_i-\bar{x}\right)^2/(n-1)}
,
and \rho=\sigma^2_{\beta}/\left(\sigma^2_{\beta}+\sigma^2\right)
is the intra-class correlation.
Value
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
.
Author(s)
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>
References
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
Examples
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)