Modified from: R CLUB 5/15/14: Power Analyses in R

Original creator: Brendan Ostlund

REFRESHER: POWER is the probability of correctly rejecting a null hypothesis that is actually false.

P(rejecting NULL hypothesis|NULL hypothesis is false)

In other words, its the prob. of finding an effect that is really present.

Power = 1 – P(Type II error)

Power, effect size, sample size, and the significance level are inter-related, and if you know 3 of these quantities you can calculate the fourth (exciting eh?). Let’s get to it.

#First, we'll need to install the "pwr" package.

# install.packages("pwr")

library(pwr)

## Warning: package 'pwr' was built under R version 3.2.2

#POWER ANALYSIS FUNCTIONS!

# T-TESTS!

Function:

pwr.t.test(n = , d = , sig.level = , power = , type = c(“two.sample”, “one.sample”, “paired”), alternative= c(“two.sided”, “less”, “greater”))

Remember, power, effect size, sample size, and the significance level are related so we can play around with the information we have to determine the missing quantity.

#For example, what is the power for a two-sided ind. samples t-test given n=25 and a medium effect size?

pwr.t.test(n=25, d=.5, sig.level= .05)

## 
##      Two-sample t test power calculation 
## 
##               n = 25
##               d = 0.5
##       sig.level = 0.05
##           power = 0.4101003
##     alternative = two.sided
## 
## NOTE: n is number in *each* group

#Another example: What is the sample size needed to achieve a pwr=.80 and a small effect size?

pwr.t.test(d=.2, sig.level= .05, power=.8)

## 
##      Two-sample t test power calculation 
## 
##               n = 393.4057
##               d = 0.2
##       sig.level = 0.05
##           power = 0.8
##     alternative = two.sided
## 
## NOTE: n is number in *each* group

# Note: if you add "$n" to the end of the command only the necessary sample size will be displayed.

pwr.t.test(d=.2, sig.level= .05, power=.8)$n

## [1] 393.4057

# R can also calculate pwr for two samples with unequal sample sizes:

pwr.t2n.test(n1= 35, n2= 17, d=.8, sig.level= .05, alternative = "greater")

## 
##      t test power calculation 
## 
##              n1 = 35
##              n2 = 17
##               d = 0.8
##       sig.level = 0.05
##           power = 0.8471659
##     alternative = greater

# ANOVA!

# what's Cohen's rules of thumb for anova effect size (i.e., "f")?

cohen.ES(test = c("anov"), size = c("small"))

## 
##      Conventional effect size from Cohen (1982) 
## 
##            test = anov
##            size = small
##     effect.size = 0.1

cohen.ES(test = c("anov"), size = c("medium"))

## 
##      Conventional effect size from Cohen (1982) 
## 
##            test = anov
##            size = medium
##     effect.size = 0.25

cohen.ES(test = c("anov"), size = c("large"))

## 
##      Conventional effect size from Cohen (1982) 
## 
##            test = anov
##            size = large
##     effect.size = 0.4

One-way ANOVA pwr analysis function:

pwr.anova.test(k = , n = , f = , sig.level = , power = )

…where k is the number of groups and n is the number of obs. IN EACH GROUP.

pwr.anova.test(k =4, n=32, f=.3, sig.level =.05)

## 
##      Balanced one-way analysis of variance power calculation 
## 
##               k = 4
##               n = 32
##               f = 0.3
##       sig.level = 0.05
##           power = 0.8101954
## 
## NOTE: n is number in each group

CORRELATION!

Function:

pwr.r.test(n = , r = , sig.level = , power = )

CHI-SQUARE

Function:

pwr.chisq.test(w =, N = , df = , sig.level =, power = )

…where “w”= effect size and “n”=number of observations.

# Quick example:

pwr.chisq.test(w=.3, N=120, df=9, sig.level=.05)

## 
##      Chi squared power calculation 
## 
##               w = 0.3
##               N = 120
##              df = 9
##       sig.level = 0.05
##           power = 0.6045631
## 
## NOTE: N is the number of observations

REGRESSION!

Function:

pwr.f2.test(u =, v = , f2 = , sig.level = , power = )

where “u”= numerator degrees of freedom (number of continuous variables + number of dummy codes – 1)

…and “v”=denominator (error) degrees of freedom

For example, find the power for a multiple regression test with 2 continuous predictors and 1 categorical

# predictor (i.e. Marital status with k=3 so 3-1=2 dummy codes) that has a large effect size and a sample size of 30.

pwr.f2.test(u =3, v =30, f2 =.35, sig.level =.05)

## 
##      Multiple regression power calculation 
## 
##               u = 3
##               v = 30
##              f2 = 0.35
##       sig.level = 0.05
##           power = 0.7807726

# Generating a table of sample sizes

# What is the required sample size for an independent samples t-test with pwr=.80?

seq <- seq(.1,.8, .1)

FindN <- rep(0,8)

for(i in 1:8) FindN[i]=pwr.t.test(d=seq[i],power=.8,sig.level=.05,type="two.sample")$n

data.frame(d=seq , N=ceiling(FindN))

##     d    N
## 1 0.1 1571
## 2 0.2  394
## 3 0.3  176
## 4 0.4  100
## 5 0.5   64
## 6 0.6   45
## 7 0.7   34
## 8 0.8   26

# can swap in other functions, just make sure arguments are relevant to test

seq <- seq(.1,.8, .1)

FindN <- rep(0,8)

for(i in 1:8) FindN[i]=pwr.anova.test(k=4, f=seq[i], power=.8,sig.level=.05)$n

data.frame(f=seq , N=ceiling(FindN))

##     f   N
## 1 0.1 274
## 2 0.2  70
## 3 0.3  32
## 4 0.4  19
## 5 0.5  12
## 6 0.6   9
## 7 0.7   7
## 8 0.8   6

# NOTE: "ceiling" is an argument for rounding numbers so the values are slightly different than if we

# computed the sample size based on specified values. For example...

pwr.anova.test(k =4,f=.3, sig.level =.05, power=.8)

## 
##      Balanced one-way analysis of variance power calculation 
## 
##               k = 4
##               n = 31.27917
##               f = 0.3
##       sig.level = 0.05
##           power = 0.8
## 
## NOTE: n is number in each group

#HELPFUL RESOURCE:

# Quick-R blog: http://www.statmethods.net/stats/power.html