1.3 Interactions and non-linear effects

1.3.1 Interactions

\[ y_i = \beta_0 + \beta_1 x_{1i} + \beta_2 x_{2i} + \beta_3 x_{1i}x_{2i} + \epsilon_i \]

We can specify the interaction between two predictors by adding an interactions list to the parameters list. Interactions can then be specified between two named variables using “:”. Interactions can be specified between two predictors at the same level, or at different hierarchical levels.

squid_data <- simulate_population(
  n=2000,
  response_name = "body_mass",
  parameters=list(
    observation=list(
      names=c("temperature","rainfall"),
      beta = c(0.5,0.3)
    ),
    residual=list(
      vcov=0.3
    ),
    interactions=list(
      names=c("temperature:rainfall"),
      beta = c(-0.1)
    )
  )
)

data <- get_population_data(squid_data)
head(data)

##     body_mass temperature     rainfall   residual temperature:rainfall
## 1  0.10157809 -0.52578806 -0.900467583 0.68195791          0.473455102
## 2  1.27037941  0.15685679  0.008508931 1.18953180          0.001334684
## 3  1.09260136  1.43214672 -0.599461006 0.47051469         -0.858516115
## 4  0.88371868  0.95614001 -1.155390758 0.64179436         -1.104715333
## 5 -0.04515612  0.05988423 -0.439114402 0.05400648         -0.026296028
## 6  2.11519536  2.73482346 -0.855646194 0.77047337         -2.340041281
##   squid_pop
## 1         1
## 2         1
## 3         1
## 4         1
## 5         1
## 6         1

coef(lm(body_mass ~ temperature * rainfall, data))

##          (Intercept)          temperature             rainfall 
##         -0.003407952          0.494634096          0.321720430 
## temperature:rainfall 
##         -0.101798139

1.3.2 Non-linear effects

Polynomial (quadratic, cubic, etc) functions are essentially interactions with the same predictor. They can therefore be specified in the same way:

squid_data <- simulate_population(
  n=2000,
  response_name = "body_mass",
  parameters=list(
    observation=list(
      names=c("temperature"),
      beta = c(0.5)
    ),
    interactions=list(
      names=c("temperature:temperature"),
      beta = c(-0.3)
    ),
    residual=list(
      vcov=0.3
    )
  )
)
data <- get_population_data(squid_data)

plot(body_mass ~ temperature, data, pch=19, cex=0.5, col=alpha(1,0.5))

coef(lm(body_mass ~ temperature + I(temperature^2), data))

##      (Intercept)      temperature I(temperature^2) 
##       0.01706406       0.49566769      -0.29948675