Douglas Curran-Everett
Explorations in statistics: power
Advan Physiol Educ 2010 34:41-43;
doi:10.1152/advan.00001.2010
#
# -- Advances_Statistics_Code_Power.R: R code for Explorations in Statistics ---
# -- Define population parameters, sample numbers, and plotting ranges ---------
#
PopMean_H0 <- 0 # population mean if null hypothesis H0 is true
PopSD <- 1 # population standard deviation
Effect_Size <- PopSD / 2 # effect size want to be able to detect
PopMean_H1 <- PopMean_H0 - Effect_Size # population mean if alternative hypothesis H1 is true
nObs <- 9 # number of observations in the sample
Alpha <- 0.10 # error rate (critical significance level)
Conf_Level <- 1 - Alpha # confidence level
# -- The Simulations -----------------------------------------------------------
#
# TheData <- round ( rnorm ( nObs, mean = PopMean_H1, sd = PopSD ), 3 ) # observations from population if H1 true
TheData <- c ( 1.894, -1.694, -0.001, -1.726, -0.727, -1.859, -0.326, -1.479, 0.476 )
TheData
## [1] 1.894 -1.694 -0.001 -1.726 -0.727 -1.859 -0.326 -1.479 0.476
Sample_Mean <- round ( mean ( TheData ), 3 )
print ( Sample_Mean )
## [1] -0.605
Sample_SD <- round ( sd ( TheData ), 3 )
print ( Sample_SD )
## [1] 1.257
t.test ( TheData, # compare observations to mean if H0 true
mu = PopMean_H0,
conf.level = Conf_Level
)
##
## One Sample t-test
##
## data: TheData
## t = -1.443, df = 8, p-value = 0.1871
## alternative hypothesis: true mean is not equal to 0
## 90 percent confidence interval:
## -1.3841 0.1748
## sample estimates:
## mean of x
## -0.6047
# -- The Simulations: power in mythical physiology paper -----------------------
#
power.t.test ( n = nObs, # number of observations
delta = PopMean_H1 - PopMean_H0, # true difference
sd = PopSD, # standard deviation
sig.level = Alpha, # error rate alpha
type = 'one.sample' # one- or two-sample problem
)
##
## One-sample t test power calculation
##
## n = 9
## delta = 0.5
## sd = 1
## sig.level = 0.1
## power = 0.3928
## alternative = two.sided
# -- The Simulations: power in mythical physiology paper -----------------------
#
power.t.test ( n = nObs,
delta = PopMean_H1 - PopMean_H0,
sd = PopSD,
sig.level = 0.01,
type = 'one.sample'
)
##
## One-sample t test power calculation
##
## n = 9
## delta = 0.5
## sd = 1
## sig.level = 0.01
## power = 0.08662
## alternative = two.sided
# -- The Simulations: number of subjects needed in mythical physiology paper ---
#
power.t.test ( power = 0.80,
delta = PopMean_H1 - PopMean_H0,
sd = PopSD,
sig.level = Alpha,
type = 'one.sample'
)
##
## One-sample t test power calculation
##
## n = 26.14
## delta = 0.5
## sd = 1
## sig.level = 0.1
## power = 0.8
## alternative = two.sided
# -- Table 1: impact of the error rate alpha on power --------------------------
#
power.t.test ( sig.level = c ( 0.01, 0.05, 0.10 ),
n = nObs,
delta = Effect_Size,
sd = PopSD,
type = 'one.sample'
)
##
## One-sample t test power calculation
##
## n = 9
## delta = 0.5
## sd = 1
## sig.level = 0.01, 0.05, 0.10
## power = 0.08662, 0.26227, 0.39277
## alternative = two.sided
# -- Table 1: impact of the difference delta mu on power -----------------------
#
power.t.test ( delta = c ( 0.5*Effect_Size, Effect_Size, 1.5*Effect_Size ),
n = nObs,
sd = PopSD,
sig.level = Alpha,
type = 'one.sample'
)
##
## One-sample t test power calculation
##
## n = 9
## delta = 0.25, 0.50, 0.75
## sd = 1
## sig.level = 0.1
## power = 0.1693, 0.3928, 0.6590
## alternative = two.sided
# -- Table 1: impact of the variability sigma on power -------------------------
#
power.t.test ( sd = c ( 2*PopSD, PopSD, 0.5*PopSD ),
n = nObs,
delta = Effect_Size,
sig.level = Alpha,
type = 'one.sample'
)
##
## One-sample t test power calculation
##
## n = 9
## delta = 0.5
## sd = 2.0, 1.0, 0.5
## sig.level = 0.1
## power = 0.1693, 0.3928, 0.8618
## alternative = two.sided
# -- Table 1: impact of the sample size n on power -----------------------------
#
power.t.test ( n = c ( nObs, 2*nObs, 4*nObs ),
delta = Effect_Size,
sd = PopSD,
sig.level = Alpha,
type = 'one.sample'
)
##
## One-sample t test power calculation
##
## n = 9, 18, 36
## delta = 0.5
## sd = 1
## sig.level = 0.1
## power = 0.3928, 0.6521, 0.9026
## alternative = two.sided