knitr::opts_chunk$set(echo = TRUE)

This vignette describes how to run power analyses for 2-way interactions that include covariates. It assumes you are already familiar with power analyses for 2-way interactions. If you aren't, check out our tutorial paper or the main vignette.

Introduction

A two-way interaction analysis with two covariates (for example, but any number is allowable) take the form:

$$ Y \sim \beta_0 + X_1\beta_1 + X_2\beta_2 + X_1X_2\beta_3+ C_1\beta_5 + C_2\beta_6 + C_1X_1\beta_7 + C_1X_2\beta_8 + C_2X_1\beta_9 + C_2X_2\beta_{10} + \epsilon $$

Where $C_i$ are the covariates, and $C_iX_i$ are interactions between the covariates and the main variables of interest. The inclusion of the $C_iX_i$ terms may surprise some users. However, it is well established that failure to include them can lead to false-positive results, because of omitted variable bias. This can occur not only when $C_iX_i$ is independent predictor of $Y$, but also simply when $cor(X_i,C_i)$ is non-zero (i.e., when a covariate is correlated with either of the main interaction terms). A few relevant citations: Hull et al., 1992, Yzerbyt et al., 2004, Keller, 2014, and blogs: Baranger, and Simonsohn.

Running an analysis

First, the function generate.interaction.cov.input() is used to setup the input correlations for the power analysis. The only input is the number of covariates. This generates a named list. All correlations are by default 0 and all reliabilites are 1.

library(InteractionPoweR)
power.input = generate.interaction.cov.input(c.num = 2) # number of covariates
head(power.input$correlations)
head(power.input$reliability)

Modify the list as-needed for your own power analysis:

power.input$correlations$r.y.x1 = .3
power.input$correlations$r.y.x2 = .2
power.input$correlations$r.y.c1 = .1
power.input$correlations$r.y.c2 = .2
power.input$correlations$r.y.x1x2 = seq(0.01,.1,.01) # inputs can be single values or vectors
power.input$correlations$r.x1.x2 = .2
power.input$correlations$r.c1.c2 = .3
power.input$correlations$r.x1.c1 = .1
power.input$correlations$r.x2.c1 = .1

Now we run the power analysis with the power_interaction_r2_covs() function:

pwr.analysis = power_interaction_r2_covs(cov.input =power.input, # our parameters
                          N = c(1000)
                          )
pwr.analysis

As always, we can plot the power curve:

plot_power_curve(power_data = pwr.analysis,x = "r.y.x1x2")

And we can estimate where the power curve intersects some desired level of power:

power_estimate(power_data = pwr.analysis,x = "r.y.x1x2",power_target = 0.8)

Even though analytic power analyses are fast, with so many parameters we can easily reach a point where the analysis will take a long time to run. Say for instance we are interested in testing out 10 different values of 4 different parameters - that's 10,000 power analyses. At that point, it can be helpful to run things in parallel. Parallel analyses can be run using the cl flag. cl=4 is probably a reasonable value for most personal computers.



dbaranger/InteractionPoweR documentation built on July 12, 2024, 3:08 p.m.