Covariates

IJLR: It's Just a Linear Regression...

All ANOVA models covered so far are linear regression model.

The latter says that

$$\underset{\text{average response}l}{\text{E}\left(Y_i{\beta }_p\right)}=\underset{\text{linear (i.e., additive) combination of explanatories}}{{\beta }_0+{\beta }_1{X}_{1i}+\cdots +{\beta }_p{X}_{pi}}$$

In an ANOVA, the model matrix $X$ simply includes columns with $-1$, $0$ and $1$ for group indicators that enforce sum-to-zero constraints.

What's in a model?

In experimental designs, the explanatories are

• experimental factors (categorical)
• continuous (dose-response)

ANCOVA = Analysis of covariance

• Analysis of variance with added continuous covariate(s) to reduce experimental error (similar to blocking).
• These continuous covariates are typically concomitant variables (measured alongside response).
• Including them in the mean response (as slopes) can help reduce the experimental error (residual error).

Control to gain power!

• enhance balance of design (randomization)
• reduce mean squared error of residuals to increase power

These steps should in principle increase power if the variables used as control are correlated with the response.

Example

Abstract of van Stekelenburg et al. (2021)

In three experiments with more than 1,500 U.S. adults who held false beliefs, participants first learned the value of scientific consensus and how to identify it. Subsequently, they read a news article with information about a scientific consensus opposing their beliefs. We found strong evidence that in the domain of genetically engineered food, this two-step communication strategy was more successful in correcting misperceptions than merely communicating scientific consensus.

Experiment 2: Genetically Engineered Food

We focus on a single experiment; preregistered exclusion criteria led to $n=442$ total sample size (unbalanced design).

Three experimental conditions:

Boost Boost Plus Consensus only (consensus)

Model formulation

Use post as response variable and prior beliefs as a control variable in the analysis of covariance.

their response was measured on a visual analogue scale ranging from –⁠100 (I am 100% certain this is false) to 100 (I am 100% certain this is true) with 0 (I don’t know) in the middle.

Plot of post vs prior response

Model formulation

Average for the $r$th replication of the $i$th experimental group is $$\begin{array}{cc}E\left({\text{post}}_{ir}\right)& =\mu +{\alpha }_i{\text{condition}}_i+\beta {\text{prior}}_{ir}.\\ Va\left({\text{post}}_{ir}\right)& ={\sigma }^2\end{array}$$

We assume that there is no interaction between condition and prior

• the slopes for prior are the same for each condition group.
• the effect of prior is linear

Contrasts of interest

1. Difference between average boosts (Boost and BoostPlus) and control (consensus)
2. Difference between Boost and BoostPlus (pairwise)

Inclusion of the prior score leads to increased precision for the mean (reduces variability).

Contrasts with ANCOVA

• The estimated marginal means will be based on a fixed value of the covariate (rather than detrended values)
• In the emmeans package, the average of the covariate is used as value.
• the difference between levels of condition are the same for any value of prior (parallel lines), but the uncertainty changes.

Multiple testing adjustments:

• Methods of Bonferroni (prespecified number of tests) and Scheffé (arbitrary contrasts) still apply
• Can't use Tukey anymore (adjusted means are not independent anymore).

Data analysis - loading data

library(emmeans)
options(contrasts = c("contr.sum", "contr.poly"))
data(SSVB21_S2, package = "hecedsm")
# Check balance
with(SSVB21_S2, table(condition))Copy Code

## condition
##     Boost BoostPlus consensus 
##       149       147       146Copy Code

Data analysis - scatterplot

library(ggplot2)
ggplot(data = SSVB21_S2,
       aes(x = prior, y = post)) + 
       geom_point() + 
       geom_smooth(method = "lm",
                   se = FALSE)Copy Code

Data analysis - model

# Check that the data are well randomized
car::Anova(lm(prior ~ condition, data = SSVB21_S2), type = 2)
# Fit linear model with continuous covariate
model1 <- lm(post ~ condition + prior, data = SSVB21_S2)
# Fit model without for comparison
model2 <- lm(post ~ condition, data = SSVB21_S2)
# Global test for differences
car::Anova(model1)
car::Anova(model2)Copy Code

Data analysis - ANOVA table

Data analysis - contrasts

emm1 <- emmeans(model1, specs = "condition")
# Note order: Boost, BoostPlus, consensus
emm2 <- emmeans(model2, specs = "condition")
# Not comparable: since one is detrended and the other isn't
contrast_list <- list(
   "boost vs control" = c(0.5,  0.5, -1), 
   #av. boosts vs consensus
   "Boost vs BoostPlus" = c(1, -1,  0))
contrast(emm1, 
         method = contrast_list, 
         p.adjust = "holm")Copy Code

Data analysis - t-tests

Data analysis - assumption checks

# Test equality of variance
levene <- car::leveneTest(
   resid(model1) ~ condition, 
   data = SSVB21_S2,
   center = 'mean')
# Equality of slopes (interaction)
car::Anova(lm(post ~ condition * prior, 
           data = SSVB21_S2),
           model1, type = 2)Copy Code

The kitchen sink approach

Should we control for more stuff?

NO! ANCOVA is a design to reduce error

• Randomization should ensure that there is no confounding
• No difference (on average) between group given a value of the covariate.
• If it isn't the case, adjustment matters

Equal trends

Testing equal slope

Compare two nested models

• Null $H_0$: model with covariate
• Alternative $H_a$: model with interaction covariate * experimental factor

Use anova to compare the models in R.

Moderator

A moderator $W$ modifies the direction or strength of the effect of an explanatory variable $X$ on a response $Y$ (interaction term).

Directed acyclic graph of moderation

Moderation in a linear regression model

In a regression model, we simply include an interaction term to the model between $W$ and $X$.

For example, if $X$ is categorical with $K$ levels and $W$ is binary or continuous, imposing sum-to-zero constraints for ${\alpha }_1,\dots ,{\alpha }_K$ and ${\beta }_1,\dots ,{\beta }_K$ gives $$\underset{\text{average response of group}k\text{at}w}{E(Y\mid X=k,W=w)}=\underset{\text{intercept of group}k}{{\alpha }_0+{\alpha }_k}+\underset{\text{slope of group}k}{({\beta }_0+{\beta }_k)}w$$

Testing for the interaction

Test jointly whether coefficients associated to $XW$ are zero, i.e., $${\beta }_1=\cdots ={\beta }_K=0.$$

The moderator $W$ can be continuous or categorical with $L\ge 2$ levels

The degrees of freedom (additional parameters for the interaction) in the $F$ test are

• $K-1$ for continuous $W$
  • are slopes parallel?
• $(K-1)\times (L-1)$ for categorical $W$
  • are all subgroup averages the same?

Example

We consider data from Garcia et al. (2010), a study on gender discrimination. Participants were given a fictional file where a women was turned down promotion in favour of male colleague despite her being clearly more experimented and qualified.

The authors manipulated the decision of the participant, with choices:

Model fit

We fit the linear model with the interaction.

data(GSBE10, package = "hecedsm")
lin_moder <- lm(respeval ~ protest*sexism, 
               data = GSBE10)
summary(lin_moder) # coefficients
car::Anova(lin_moder, type = 2) # testsCopy Code

ANOVA table

Effects

Comparisons between groups

Simple effects and comparisons must be done for a fixed value of sexism (since the slopes are not parallel).

The default value in emmeans is the mean value of sexism, but we could query for averages at different values of sexism (below for empirical quartiles).

quart <-  quantile(GSBE10$sexism, probs = c(0.25, 0.5, 0.75))
emmeans(lin_moder, 
        specs = "protest", 
        by = "sexism",
        at = list("sexism" = quart))Copy Code

Sensitivity analysis

The Johnson and Neyman (1936) method looks at the range of values of moderator $W$ for which difference between treatments (binary $X$) is not statistically significant.

lin_moder2 <- lm(
  respeval ~ protest*sexism, 
  data = GSBE10 |> 
  # We dichotomize the manipulation, pooling protests together
  dplyr::mutate(protest = as.integer(protest != "no protest")))
# Test for equality of slopes/intercept for two protest groups
anova(lin_moder, lin_moder2)
# p-value of 0.18: fail to reject individual = collective.Copy Code

Syntax for plot

jn <- interactions::johnson_neyman(
  model = lin_moder2, # linear model
  pred = protest, # binary experimental factor
  modx = sexism, # moderator
  control.fdr = TRUE, # control for false discovery rate
  mod.range = range(GSBE10$sexism)) # range of values for sexism
jn$plotCopy Code

Plot of Johnson−Neyman intervals

Johnson−Neyman plot for difference between protest and no protest as a function of sexism.

Moderation

More generally, moderation refers to any explanatory variable (whether continuous or categorical) which interacts with the experimental manipulation.

• For categorical-categorical, this is a multiway ANOVA model
• For continuous-categorical, use linear regression

Summary