Task 1
Use data for Germany and the variable coding for model3 and model4 from the tutorial (i.e., both for the cleaning & recode, as well as the linear model formulas). Do a modelsummary table displaying both model outputs. Interpret the coefficients for the MLR without the interaction (i.e., the first displayed model), and interpret model fit metrics for both models.
germany_data <- read_fst("germany_data.fst")
df <- germany_data
df <- df %>%
mutate(behave = ipbhprp,
secure = impsafe,
safety = ipstrgv,
tradition = imptrad,
rules = ipfrule) %>%
mutate(across(c("behave", "secure", "safety", "tradition", "rules"),
~ na_if(.x, 7) %>% na_if(8) %>% na_if(9))) %>%
# Apply the reverse coding
mutate(across(c("behave", "secure", "safety", "tradition", "rules"), ~ 7 - .x ))
df$auth <- scales::rescale(df$behave +
df$secure +
df$safety +
df$tradition +
df$rules, to=c(0,100), na.rm=TRUE)
df <- df %>%
mutate(
polID = case_when(
lrscale %in% 0:3 ~ "Left",
lrscale %in% 7:10 ~ "Right",
lrscale %in% 4:6 ~ "Moderate",
lrscale %in% c(77, 88, 99) ~ NA_character_
),
religious = case_when(
rlgdgr %in% c(77, 88, 99) ~ NA_real_,
TRUE ~ rlgdgr
)
)
df <- df %>% filter(!is.na(auth))
model3 <- lm(auth ~ polID + religious, data = df, weights = weight)
model4 <- lm(auth ~ polID + religious + polID*religious, data = df, weights = weight)
modelsummary(
list(model3, model4),
fmt = 1,
estimate = c( "{estimate} ({std.error}){stars}",
"{estimate} ({std.error}){stars}"),
statistic = NULL)
| Â (1) | Â Â (2) | |
|---|---|---|
| (Intercept) | 57.5 (0.7)*** | 59.6 (0.9)*** |
| polIDModerate | 5.3 (0.7)*** | 2.6 (1.1)* |
| polIDRight | 9.3 (1.1)*** | 3.4 (1.9)+ |
| religious | 0.3 (0.1)** | −0.3 (0.2)+ |
| polIDModerate × religious | 0.8 (0.2)** | |
| polIDRight × religious | 1.5 (0.4)*** | |
| Num.Obs. | 2855 | 2855 |
| R2 | 0.035 | 0.041 |
| R2 Adj. | 0.034 | 0.039 |
| AIC | 24680.0 | 24665.5 |
| BIC | 24709.7 | 24707.2 |
| Log.Lik. | −12334.978 | −12325.764 |
| RMSE | 17.13 | 17.13 |
Task 1 Answer: PollDModerate (5.3) suggests that people who identify as politically moderate are 5.3 units higher than the left (the reference category). PollDRight (9.3) shows that people who identify as right (politically) are 9.3 units higher than people who identify as more left-sided (on average).The religious predictor is 0.3, which means that with every single-unit increase in religiousness, the associated expected average will also increase by 0.3 (on the authoritarian values scale).The *** for each coefficient tells us that the p-value is less than 0.001, which means that the probability of the observed result being due to chance is less than 0.1%. The null hypothesis is incorrectly rejected.The fit metric, adjusted R-squared, is 0.034(less than 1) which means that the model has not an ideal predictor.
Task 2 Now generate the model4 interaction plot that we did in the tutorial, but again using the German data instead of the French. Interpret.
germany_data <- read_fst("germany_data.fst")
interaction_plot <- effect("polID*religious", model4, na.rm=TRUE)
plot(interaction_plot,
main="Interaction effect",
xlab="Religiousness",
ylab="Authoritarian attitudes scale")
interaction_plot
##
## polID*religious effect
## religious
## polID 0 2 5 8 10
## Left 59.61010 58.93398 57.91980 56.90562 56.22950
## Moderate 62.21299 63.14317 64.53843 65.93370 66.86388
## Right 63.00456 65.26702 68.66072 72.05442 74.31689
Task 2 Answer: The downward slope of pollD Left shows that according to the Authoritarian attitudes scale, more people that identify as politically left also identify as less religious. The slight upward slope of PollD Moderate shows that people who identify as politically moderate also identify as slightly more religious. Lastly, the steep upward slope of pollD Right shows us that according to the Authoritarian attitudes scale, more individuals that identify as politically right also have a stronger sense of religiosity.
Task 3 Use data for the Netherlands and the variable coding for model5 and model6 from the tutorial (i.e., both for the cleaning & recode, as well as the linear model formulas). Do a modelsummary table displaying both model outputs. Interpret the coefficients for the MLR without the interaction (i.e., the first displayed model), and interpret model fit metrics for both models.
netherlands_data <- read_fst("netherlands_data.fst")
df <- df %>%
mutate(
polID = case_when(
lrscale %in% 0:3 ~ "Left",
lrscale %in% 7:10 ~ "Right",
lrscale %in% 4:6 ~ "Moderate",
lrscale %in% c(77, 88, 99) ~ NA_character_
),
religious = case_when(
rlgdgr %in% c(77, 88, 99) ~ NA_real_,
TRUE ~ rlgdgr
)
)
df <- df %>%
mutate(
cohort = ifelse(yrbrn < 1930 | yrbrn > 2000, NA, yrbrn),
# Recoding generational cohorts based on the year of birth (yrbrn).
# The year of birth is categorized into different generational cohorts.
# Interwar (1900-1945), Baby Boomers (1946-1964), Gen X (1965-1979), Millennials (1980-1996).
# The 'TRUE' line is a catch-all that keeps the original year of birth for those not in these ranges.
gen = case_when(
yrbrn %in% 1900:1945 ~ "1",
yrbrn %in% 1946:1964 ~ "2",
yrbrn %in% 1965:1979 ~ "3",
yrbrn %in% 1980:1996 ~ "4",
TRUE ~ as.character(yrbrn)
),
# After recoding, the gen variable is converted into a factor with labels for clearer interpretation.
# Factors are used in R to handle categorical variables.
gen = factor(gen,
levels = c("1", "2", "3", "4"),
labels = c("Interwar", "Baby Boomers", "Gen X", "Millennials"))
)
table(df$gen)
##
## Interwar Baby Boomers Gen X Millennials
## 5360 8583 5601 4652
df <- df %>%
mutate(religion = case_when(
rlgblg == 2 ~ "No",
rlgblg == 1 ~ "Yes",
rlgblg %in% c(7, 8, 9) ~ NA_character_,
TRUE ~ as.character(rlgblg)
))
# check
table(df$religion)
##
## No Yes
## 11088 13750
df <- df %>%
mutate(ID = case_when(
lrscale >= 0 & lrscale <= 4 ~ "Left",
lrscale >= 6 & lrscale <= 10 ~ "Right",
lrscale > 10 ~ NA_character_, # Set values above 10 as NA
TRUE ~ NA_character_ # Ensure value 5 and any other unexpected values are set as NA
))
table(df$ID)
##
## Left Right
## 9745 4882
model5 <- lm(auth ~ religion + ID + gen, data = df, weights = weight)
model6 <- lm(auth ~ religion + ID + gen + religion*gen, data = df, weights = weight)
modelsummary(
list(model5, model6),
fmt = 1,
estimate = c( "{estimate} ({std.error}){stars}",
"{estimate} ({std.error}){stars}"),
statistic = NULL,
coef_omit = "Intercept")
| Â (1) | Â Â (2) | |
|---|---|---|
| religionYes | 1.3 (0.8) | −0.9 (2.0) |
| IDRight | 7.1 (0.9)*** | 7.3 (0.9)*** |
| genBaby Boomers | −7.6 (1.1)*** | −8.2 (1.8)*** |
| genGen X | −9.7 (1.3)*** | −11.5 (1.9)*** |
| genMillennials | −12.7 (1.3)*** | −15.7 (2.0)*** |
| religionYes × genBaby Boomers | 0.8 (2.3) | |
| religionYes × genGen X | 3.2 (2.5) | |
| religionYes × genMillennials | 5.8 (2.7)* | |
| Num.Obs. | 1761 | 1761 |
| R2 | 0.103 | 0.106 |
| R2 Adj. | 0.100 | 0.102 |
| AIC | 15122.2 | 15121.7 |
| BIC | 15160.5 | 15176.5 |
| Log.Lik. | −7554.112 | −7550.862 |
| RMSE | 16.76 | 16.76 |
interaction_plot <- effect("religion*gen", model6, na.rm=TRUE)
plot(interaction_plot,
main="Interaction effect",
xlab="Generations",
ylab="Authoritarian attitudes scale")
interaction_plot
##
## religion*gen effect
## gen
## religion Interwar Baby Boomers Gen X Millennials
## No 70.47065 62.23112 58.99075 54.73118
## Yes 69.61993 62.20566 61.30541 59.63364
cat_plot(model5, pred = gen, modx = religion, jnplot = TRUE)
## Warning: gen and religion are not included in an interaction with one another in the
## model.
cat_plot(model6, pred = gen, modx = religion, jnplot = TRUE)
Task 3 Answer:
Model 5: religionYes (1.3) suggests that people who identify as religious are 1.3 units higher than average. IDRight (7.1) shows that people who identify as politically right are 7.1 units higher than the average, and are also more religious. The adjusted R-squared is 0.100, which means that the model does not explain most of the variation in the data. The *** for each coefficient (genBaby Boomers, genGen X, genMillenials) means the p-value is less than 0.001, which means the explanatory variables are not statistically significant and the null hypothesis is incorrectly rejected. This means that the model is a poor predictor.
Model 6: religionYes (-0.9) shows us that people who identify as religious are 0.9 units below the average. IDRight (7.3) shows us that people who identify as politically right are 7.3 units above the average. the adjusted R-squared is 0.102, so the model does not account for most of the variation of the data. The *** on each coefficient shows that the p-value is 0.001, meaning they are not statistically significant.
Task 4 Produce the model7 interaction plot from the tutorial, but again using the Netherlands data. Interpret.
netherlands_data <- read_fst("netherlands_data.fst")
interaction_plot
##
## religion*gen effect
## gen
## religion Interwar Baby Boomers Gen X Millennials
## No 70.47065 62.23112 58.99075 54.73118
## Yes 69.61993 62.20566 61.30541 59.63364
model7 <- lm(auth ~ + cohort + polID + cohort*polID, data = df, weights = weight)
interact_plot(model7, pred = cohort, modx = polID, jnplot = TRUE)
Task 4 Answer: None of the lines converge, indicating that no
interaction occurs.
Task 5 Produce the model8 interaction plot from the tutorial, but again using the Netherlands data. Interpret.
model8 <- lm(auth ~ religious + ID + religious*ID, data = df, weights = weight)
interaction_plot
##
## religion*gen effect
## gen
## religion Interwar Baby Boomers Gen X Millennials
## No 70.47065 62.23112 58.99075 54.73118
## Yes 69.61993 62.20566 61.30541 59.63364
interact_plot(model8, pred = religious, modx = ID, jnplot = TRUE)
Task 5 Answer: The lines are not parallel, so there is an interaction
between identifying as politically left/right and religiosity.
Individuals who identify as politically right have stronger religious
attitudes and individuals who identify as politically left have
significantly weaker religious attitudes.