This vignette replicates several of the examples found in the G*Power manual (version 3.1). It is not meant to be exhaustive but instead demonstrates how power analyses could be computed and extended using simulation methodology by either editing the default functions found within the package or by creating a new user-defined function for yet-to-be-defined statistical analysis contexts.
Power associated with the hypotheses \[H_0:\, \rho-\rho_0=0\] \[H_1:\, \rho-\rho_0\ne 0\]
where \(\rho\) is the population correlation and \(\rho_0\) the null hypothesis constant.
Sample size estimate to reject \(H_0:\, \rho_0=.60\) in correlation analysis with \(1-\beta=.95\) probability when \(\rho=.65\).
##
## Design conditions:
##
## # A tibble: 1 × 5
## n r rho sig.level power
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 NA 0.65 0.6 0.05 0.95
##
## Estimate of n: 1931.4
## 95% Prediction Interval: [1901.1, 1958.2]
G*power estimates \(n\) to be 1929 using the same Fisher z-transformation approximation.
Unlike the previous section, this is the more canonical version of the hypotheses involving correlation coefficients. Power associated with \(\rho = .3\) with 100 pairs of observations, tested against \(\rho_0=0\).
##
## Design conditions:
##
## # A tibble: 1 × 4
## n r sig.level power
## <dbl> <dbl> <dbl> <lgl>
## 1 100 0.3 0.05 NA
##
## Estimate of power: 0.861
## 95% Confidence Interval: [0.854, 0.867]
Sample size estimate to reject \(H_0:\, \rho_0=0\) in correlation analysis with \(1-\beta=.95\) probability when \(\rho=.3\).
##
## Design conditions:
##
## # A tibble: 1 × 4
## n r sig.level power
## <dbl> <dbl> <dbl> <dbl>
## 1 NA 0.3 0.05 0.95
##
## Estimate of n: 138.3
## 95% Prediction Interval: [136.1, 140.4]
Compare to approximate result from pwr
package.
##
## approximate correlation power calculation (arctangh transformation)
##
## n = 137.8
## r = 0.3
## sig.level = 0.05
## power = 0.95
## alternative = two.sided
##
## Design conditions:
##
## # A tibble: 1 × 4
## r.ab1 r.ab2 sig.level power
## <dbl> <dbl> <dbl> <lgl>
## 1 0.75 0.88 0.05 NA
##
## Estimate of power: 0.727
## 95% Confidence Interval: [0.718, 0.735]
G*power gives power of .726.
# solution per group
out <- p_t.test(r = .25, n = NA, two.tailed=FALSE) |>
Spower(power = .95, interval=c(50, 200))
out
##
## Design conditions:
##
## # A tibble: 1 × 5
## n r two.tailed sig.level power
## <dbl> <dbl> <lgl> <dbl> <dbl>
## 1 NA 0.25 FALSE 0.05 0.95
##
## Estimate of n: 81.1
## 95% Prediction Interval: [79.4, 82.7]
## [1] 164
G*power gives the result \(N=164\).
Relatedly, one can specify \(d\), Cohen’s standardized mean-difference effect size, instead of \(r\) since \(d\) is easily converted to \(r\).
## [1] 0.4183 0.5817
## [1] 0.3978 0.6022
## [1] 0.2063 0.2589
G*power gives \(n=463\), though uses the SE value at the null (Score test). p_r.cat()
, on
the other hand, defaults to the Wald approach where the SE is at the MLE. To switch, use
score=FALSE
, though note that this requires twice as many computations.
p_r.cat(n=NA, r=0.2399846, tauX=tauX, tauY=tauY,
score=FALSE, two.tailed=FALSE) |>
Spower(power = .95, interval=c(100, 500),
parallel=TRUE, check.interval=FALSE)
##
## Design conditions:
##
## # A tibble: 1 × 8
## n r tauX tauY score two.tailed sig.level power
## <dbl> <dbl> <dbl> <dbl> <lgl> <lgl> <dbl> <dbl>
## 1 NA 0.240 0.206 0.259 FALSE FALSE 0.05 0.95
##
## Estimate of n: 462.9
## 95% Prediction Interval: [458.5, 466.6]
One sample, one-tailed proportions test given data generated from a population tested against \(\pi_0 = .65\) with \(g=.15\) (hence, \(\pi = .65+g=.80\)) with \(n=20\).
##
## Design conditions:
##
## # A tibble: 1 × 4
## n two.tailed sig.level power
## <dbl> <lgl> <dbl> <lgl>
## 1 20 FALSE 0.05 NA
##
## Estimate of power: 0.423
## 95% Confidence Interval: [0.413, 0.432]
G*power gives the estimate \(1-\beta=.4112\).
Fisher’s exact test is also supported by using the argument exact = TRUE
.
##
## Design conditions:
##
## # A tibble: 1 × 5
## n two.tailed exact sig.level power
## <dbl> <lgl> <lgl> <dbl> <lgl>
## 1 20 FALSE TRUE 0.05 NA
##
## Estimate of power: 0.415
## 95% Confidence Interval: [0.405, 0.425]
Normal as the parent distribution.
##
## Design conditions:
##
## # A tibble: 1 × 6
## n d type two.tailed sig.level power
## <dbl> <dbl> <chr> <lgl> <dbl> <lgl>
## 1 649 0.1 one.sample FALSE 0.05 NA
##
## Estimate of power: 0.794
## 95% Confidence Interval: [0.786, 0.802]
G*power gives the power estimate of .800.
A one-sample sign test with Laplace distribution as the parent:
library(VGAM)
# generate data with scale 0-1 for d effect size to be same as mean
parent <- function(n, d, scale=sqrt(1/2))
VGAM::rlaplace(n, d, scale=scale)
p_wilcox.test(n=11, d=.8, parent1=parent, type='one.sample',
two.tailed=FALSE, correct = FALSE) |> Spower()
##
## Design conditions:
##
## # A tibble: 1 × 7
## n d type correct two.tailed sig.level power
## <dbl> <dbl> <chr> <lgl> <lgl> <dbl> <lgl>
## 1 11 0.8 one.sample FALSE FALSE 0.05 NA
##
## Estimate of power: 0.809
## 95% Confidence Interval: [0.801, 0.817]
G*power gives the estimate .830.
obrien2002 <- matrix(c(.54, .32, .08, .06), 2, 2,
dimnames = list('Treatment' = c('Yes', 'No'),
'Standard' = c('Yes', 'No')))
obrien2002
## Standard
## Treatment Yes No
## Yes 0.54 0.08
## No 0.32 0.06
##
## Design conditions:
##
## # A tibble: 1 × 4
## n two.tailed sig.level power
## <dbl> <lgl> <dbl> <lgl>
## 1 50 FALSE 0.05 NA
##
## Estimate of power: 0.838
## 95% Confidence Interval: [0.834, 0.842]
G*Power gives .839 (\(\alpha = .032\)). Slightly more power can be achieved when not using the continuity correction, though in general this is not recommended in practice.
##
## Design conditions:
##
## # A tibble: 1 × 5
## n two.tailed correct sig.level power
## <dbl> <lgl> <lgl> <dbl> <lgl>
## 1 50 FALSE FALSE 0.05 NA
##
## Estimate of power: 0.887
## 95% Confidence Interval: [0.881, 0.893]
##
## Design conditions:
##
## # A tibble: 1 × 6
## n R2 k fixed sig.level power
## <dbl> <dbl> <dbl> <lgl> <dbl> <lgl>
## 1 95 0.1 5 FALSE 0.05 NA
##
## Estimate of power: 0.667
## 95% Confidence Interval: [0.658, 0.676]
G*power gives 0.662 using a one-tailed test criterion.
##
## Design conditions:
##
## # A tibble: 1 × 5
## n R2 k sig.level power
## <dbl> <dbl> <dbl> <dbl> <lgl>
## 1 95 0.1 5 0.05 NA
##
## Estimate of power: 0.660
## 95% Confidence Interval: [0.651, 0.669]
G*power gives \(1-\beta = .673\).
Note that k
is total IVs, k.R2_0
is number of IVs for baseline model.
##
## Design conditions:
##
## # A tibble: 1 × 7
## n R2 k R2_0 k.R2_0 sig.level power
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <lgl>
## 1 90 0.3 9 0.25 5 0.01 NA
##
## Estimate of power: 0.231
## 95% Confidence Interval: [0.223, 0.239]
G*power gives \(1-\beta = .241\). Solving the sample size to achieve 80% power
p_lm.R2(n=NA, R2=.3, R2_0 = .25, k=9, k.R2_0=5) |>
Spower(sig.level=.01, power=.8, interval=c(100, 400))
##
## Design conditions:
##
## # A tibble: 1 × 7
## n R2 k R2_0 k.R2_0 sig.level power
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 NA 0.3 9 0.25 5 0.01 0.8
##
## Estimate of n: 243.4
## 95% Prediction Interval: [242.1, 244.7]
G*power gives \(n = 242\).
Compare model with 12 IVs to model with 9 IVs.
##
## Design conditions:
##
## # A tibble: 1 × 8
## n R2 k R2_0 k.R2_0 R2.resid sig.level power
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <lgl>
## 1 200 0.16 12 0.1 9 0.8 0.01 NA
##
## Estimate of power: 0.759
## 95% Confidence Interval: [0.751, 0.767]
G*power gives \(1-\beta = .767\).
##
## Design conditions:
##
## # A tibble: 1 × 5
## n k f sig.level power
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 NA 10 0.25 0.05 0.95
##
## Estimate of n: 38.4
## 95% Prediction Interval: [37.7, 39.1]
G*power gives the estimate \(n=39\).
Fixing \(n=200\) in total (hence, \(n=200/k=20\)) and performing a compromise analysis assuming \(q=\frac{\beta}{\alpha}=1\),
##
## Design conditions:
##
## # A tibble: 1 × 6
## n k f sig.level power beta_alpha
## <dbl> <dbl> <dbl> <dbl> <lgl> <dbl>
## 1 20 10 0.25 NA NA 1
##
## Estimate of Type I error rate (alpha/sig.level): 0.159
## 95% Confidence Interval: [0.154, 0.163]
##
## Estimate of power (1-beta): 0.841
## 95% Confidence Interval: [0.837, 0.846]
G*Power gives \(\alpha=\beta=0.159\).
Test coefficients across distinct datasets with similar form. In this case
\[Y_1 = \beta_0 + \beta_1 X_1 + \epsilon\] \[Y_2 = \beta_0 + \beta_2 X_2 + \epsilon\]
where the null of interest is
\[H_0:\, \beta_1 - \beta_2 = 0\]
To do this a multiple linear regression model is setup with three variables
\[Y = \beta_0 + \beta_1 X + \beta_2 S + \beta_3 S\cdot X + \epsilon\] where \(Y=[Y_1, Y_2]\), \(X = [X_1, X_2]\), and \(S\) is a binary indicator variable indicating whether the observations were in the second sample. Hence, the same null hypothesis can be evaluated using this augmented model testing the null \[H_0:\, \beta_3 = 0\]
We start by defining the population generating model to replace the gen_glm()
function that is the default in p_glm()
. This generating function is organized such that a data.frame
is returned with the columns y
, X
, and S
, where the interaction effect reflects the magnitude of the difference between the \(\beta\) coefficients across the independent samples.
gen_twogroup <- function(n, dbeta, sdx1, sdx2, sigma, n2_n1 = 1, ...){
X1 <- rnorm(n, sd=sdx1)
X2 <- rnorm(n*n2_n1, sd=sdx2)
X <- c(X1, X2)
N <- length(X)
S <- c(rep(0, n), rep(1, N-n))
y <- dbeta * X*S + rnorm(N, sd=sigma)
dat <- data.frame(y, X, S)
dat
}
To demonstrate, the post-hoc power for the described example in G*Power is the following.
p_glm(formula=y~X*S, test="X:S = 0",
n=28, n2_n1=44/28, sdx1=9.02914, sdx2=11.86779, dbeta=0.01592,
sigma=0.5578413, gen_fun=gen_twogroup) |> Spower()
##
## Design conditions:
##
## # A tibble: 1 × 8
## test sigma n sdx1 sdx2 dbeta sig.level power
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <lgl>
## 1 X:S = 0 0.55784 28 9.0291 11.868 0.01592 0.05 NA
##
## Estimate of power: 0.206
## 95% Confidence Interval: [0.198, 0.214]
For the a priori power analysis to achieve a power of .80
p_glm(formula=y~X*S, test="X:S = 0",
n=NA, n2_n1=44/28, sdx1=9.02914, sdx2=11.86779, dbeta=0.01592,
sigma=0.5578413, gen_fun=gen_twogroup) |>
Spower(power=.8, interval=c(100, 1000))
##
## Design conditions:
##
## # A tibble: 1 × 8
## test sigma n sdx1 sdx2 dbeta sig.level power
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 X:S = 0 0.55784 NA 9.0291 11.868 0.01592 0.05 0.8
##
## Estimate of n: 165.5
## 95% Prediction Interval: [164.1, 166.8]
G*Power gives the estimate for \(n\) to be 163 (and therefore 256 in the second group given \(n2_n1\)).
# solve n for variance ratio of 1/1.5 = 2/3, one.tailed, 80% power
p_var.test(n=NA, vars=1, sigma2=1.5, two.tailed=FALSE) |>
Spower(power=.80, interval=c(10, 200))
##
## Design conditions:
##
## # A tibble: 1 × 6
## n vars sigma2 two.tailed sig.level power
## <dbl> <dbl> <dbl> <lgl> <dbl> <dbl>
## 1 NA 1 1.5 FALSE 0.05 0.8
##
## Estimate of n: 81.1
## 95% Prediction Interval: [80.4, 81.8]
G*power gives sample size of 81.
For a two-sample equality of variance test with equal sample sizes,
# solve n for variance ratio of 1/1.5 = 2/3, two.tailed, 80% power
p_var.test(n=NA, vars=c(1, 1.5), two.tailed=TRUE) |>
Spower(power=.80, interval=c(50, 300))
##
## Design conditions:
##
## # A tibble: 1 × 4
## n two.tailed sig.level power
## <dbl> <lgl> <dbl> <dbl>
## 1 NA TRUE 0.05 0.8
##
## Estimate of n: 193.2
## 95% Prediction Interval: [190.9, 195.3]
G*Power gives estimate of 193 per group.
Estimate sample size (\(n\)) per group in independent samples \(t\)-test, one-tailed, medium effect size (\(d=0.5\)), \(\alpha=0.05\), 95% power (\(1-\beta = 0.95\)), equal sample sizes (\(\frac{n_2}{n_1}=1\)).
##
## Design conditions:
##
## # A tibble: 1 × 5
## n d two.tailed sig.level power
## <dbl> <dbl> <lgl> <dbl> <dbl>
## 1 NA 0.5 FALSE 0.05 0.95
##
## Estimate of n: 88.7
## 95% Prediction Interval: [85.7, 91.6]
G*power estimate is 88 per group, Spower
estimate is 88.676 with the 95% CI [85.7064, 91.553].
##
## Design conditions:
##
## # A tibble: 1 × 5
## n d type sig.level power
## <dbl> <dbl> <chr> <dbl> <lgl>
## 1 50 0.42164 paired 0.05 NA
##
## Estimate of power: 0.841
## 95% Confidence Interval: [0.834, 0.848]
G*power gives power estimate of .832.
p_t.test(n=NA, d=.625, two.tailed=FALSE, type='one.sample') |>
Spower(power=.95, interval=c(10, 100))
##
## Design conditions:
##
## # A tibble: 1 × 6
## n d type two.tailed sig.level power
## <dbl> <dbl> <chr> <lgl> <dbl> <dbl>
## 1 NA 0.625 one.sample FALSE 0.05 0.95
##
## Estimate of n: 29.1
## 95% Prediction Interval: [28.6, 29.5]
G*power gives sample size of \(n=30\).
##
## Design conditions:
##
## # A tibble: 1 × 5
## n d type sig.level power
## <dbl> <dbl> <chr> <dbl> <dbl>
## 1 NA 0.1 one.sample 0.01 0.9
##
## Estimate of n: 1503.4
## 95% Prediction Interval: [1490.4, 1516.4]
G*power gives sample size of \(n=1492\).
##
## Design conditions:
##
## # A tibble: 1 × 6
## n d type two.tailed sig.level power
## <dbl> <dbl> <chr> <lgl> <dbl> <lgl>
## 1 649 0.1 one.sample FALSE 0.05 NA
##
## Estimate of power: 0.795
## 95% Confidence Interval: [0.787, 0.803]
G*power provide power estimate of .800.
library(VGAM)
parent1 <- function(n, d) VGAM::rlaplace(n, d, scale=sqrt(1/2))
parent2 <- function(n, d) VGAM::rlaplace(n, scale=sqrt(1/2))
nr <- 134/67
p_wilcox.test(n=67, n2_n1=nr, d=0.375, parent1=parent1, parent2=parent2) |>
Spower()
##
## Design conditions:
##
## # A tibble: 1 × 4
## n d sig.level power
## <dbl> <dbl> <dbl> <lgl>
## 1 67 0.375 0.05 NA
##
## Estimate of power: 0.845
## 95% Confidence Interval: [0.838, 0.852]
G*power gives power of .847.