1

Let’s re-estimate our model from the previous exercise, just with another dependent variable. Again, consult the codebook that lives in the ./data folder or the full Allbus 2021 codebook and choose another variable that may make sense to be inserted into a regression.

If you’re having a hard time finding suitable variables, what about one of the dependent variables from the first set of exercises for this session:

  1. mi05 (accepting refugees from war-torn countries)
  2. ls01 (satisfaction with general living situation)
  3. fr12 (man taking care of the household and children)
Be aware that you may have to conduct some recoding.
# We chose ls01 for the following tasks/solutions.

2

Re-run your analysis just by switching the dependent variables.
library(dplyr)
library(haven)
library(parameters)

allbus_2021_dvII_3 <- 
  read_spss("./data/allbus_2021/ZA5280_v1-0-0.sav") %>% 
  mutate(life_satisfaction = as.numeric(ls01))
## Converting atomic to factors. Please wait...
linear_model_2 <-
  lm(
    life_satisfaction  ~ age + eastwest,
    data = allbus_2021_dvII_3
  )

model_parameters(linear_model_2)
## Parameter    | Coefficient |       SE |         95% CI | t(5107) |      p
## -------------------------------------------------------------------------
## (Intercept)  |        7.85 |     0.09 | [ 7.69,  8.02] |   92.19 | < .001
## age          |    9.21e-03 | 1.51e-03 | [ 0.01,  0.01] |    6.09 | < .001
## eastwest [2] |       -0.25 |     0.06 | [-0.36, -0.14] |   -4.48 | < .001
## 
## Uncertainty intervals (equal-tailed) and p-values (two-tailed) computed using a Wald t-distribution approximation.

3

Now extract the prediction data for your main independent variable for this model.
You can use the function get_model_data() from the sjPlot package. You should set the option type = "pred" and provide the name of your independent variable in the terms option.
library(sjPlot)

predictions_model_2 <-
  get_model_data(
    linear_model_2,
    type = "pred",
    terms = "age"
  )

4

Repeat the previous step for your original model and combine both datasets. The resulting data should be a data.frame or a tibble.
  • For combining the data, you simply have append both data sets rowwise. You can either use the bind_rows() function from the dplyr package or rbind() from base R.
  • Make sure that you add an indicator variable for the model for each of the data sets before during the combination.
  • You can convert any (well, a lot of…) data objects to a tibble using as_tibble() from the tibble package.
library(dplyr)
library(tibble)

linear_model<-
  lm(
    political_orientation ~ age + eastwest,
    data = allbus_2021_dvII_2
  )

predictions_model <-
  get_model_data(
    linear_model,
    type = "pred",
    terms = "age"
  )

predictions <-
  bind_rows(
    predictions_model %>% 
      mutate(model = "Model 1"),
    predictions_model_2 %>% 
      mutate(model = "Model 2")
  ) %>% 
  as_tibble()

5

Create a faceted ggplot with the regression predictions (without further customizing anything).
You can use facet_wrap() here.
library(ggplot2)

ggplot(
  predictions,
  aes(x, predicted)
  ) + 
  geom_line() +
  geom_line(aes(x, conf.low), linetype = "dashed") +
  geom_line(aes(x, conf.high), linetype = "dashed") +
  facet_wrap(~model)