Tutorial 7 (SOC252)
Sebastien Parker
2023-11-02
Setting up environment
setwd("C:/Users/Owner/Downloads/All Downloads/Stepping 2 Workspace")
getwd()## [1] "C:/Users/Owner/Downloads/All Downloads/Stepping 2 Workspace"# List of packages
# List of packages
packages <- c("tidyverse", "infer", "fst", "modelsummary", "effects", "ggplot2", "tidyverse", "modelsummary", "forcats", "RColorBrewer", "ggplot2",
"fst", "viridis", "knitr", "kableExtra", "rmarkdown", "ggridges", "viridis", "questionr") # add any you need here
# Install packages if they aren't installed already
new_packages <- packages[!(packages %in% installed.packages()[,"Package"])]
if(length(new_packages)) install.packages(new_packages)
# Load the packages
lapply(packages, library, character.only = TRUE)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.3 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.3 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.0
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
## Loading required package: carData
##
## lattice theme set by effectsTheme()
## See ?effectsTheme for details.
##
## Loading required package: viridisLite
##
##
## Attaching package: 'kableExtra'
##
##
## The following object is masked from 'package:dplyr':
##
## group_rows## [[1]]
## [1] "lubridate" "forcats" "stringr" "dplyr" "purrr" "readr"
## [7] "tidyr" "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [13] "grDevices" "utils" "datasets" "methods" "base"
##
## [[2]]
## [1] "infer" "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [7] "readr" "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [13] "graphics" "grDevices" "utils" "datasets" "methods" "base"
##
## [[3]]
## [1] "fst" "infer" "lubridate" "forcats" "stringr" "dplyr"
## [7] "purrr" "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [13] "stats" "graphics" "grDevices" "utils" "datasets" "methods"
## [19] "base"
##
## [[4]]
## [1] "modelsummary" "fst" "infer" "lubridate" "forcats"
## [6] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [11] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [16] "grDevices" "utils" "datasets" "methods" "base"
##
## [[5]]
## [1] "effects" "carData" "modelsummary" "fst" "infer"
## [6] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [11] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [16] "stats" "graphics" "grDevices" "utils" "datasets"
## [21] "methods" "base"
##
## [[6]]
## [1] "effects" "carData" "modelsummary" "fst" "infer"
## [6] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [11] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [16] "stats" "graphics" "grDevices" "utils" "datasets"
## [21] "methods" "base"
##
## [[7]]
## [1] "effects" "carData" "modelsummary" "fst" "infer"
## [6] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [11] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [16] "stats" "graphics" "grDevices" "utils" "datasets"
## [21] "methods" "base"
##
## [[8]]
## [1] "effects" "carData" "modelsummary" "fst" "infer"
## [6] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [11] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [16] "stats" "graphics" "grDevices" "utils" "datasets"
## [21] "methods" "base"
##
## [[9]]
## [1] "effects" "carData" "modelsummary" "fst" "infer"
## [6] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [11] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [16] "stats" "graphics" "grDevices" "utils" "datasets"
## [21] "methods" "base"
##
## [[10]]
## [1] "RColorBrewer" "effects" "carData" "modelsummary" "fst"
## [6] "infer" "lubridate" "forcats" "stringr" "dplyr"
## [11] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [16] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [21] "datasets" "methods" "base"
##
## [[11]]
## [1] "RColorBrewer" "effects" "carData" "modelsummary" "fst"
## [6] "infer" "lubridate" "forcats" "stringr" "dplyr"
## [11] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [16] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [21] "datasets" "methods" "base"
##
## [[12]]
## [1] "RColorBrewer" "effects" "carData" "modelsummary" "fst"
## [6] "infer" "lubridate" "forcats" "stringr" "dplyr"
## [11] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [16] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [21] "datasets" "methods" "base"
##
## [[13]]
## [1] "viridis" "viridisLite" "RColorBrewer" "effects" "carData"
## [6] "modelsummary" "fst" "infer" "lubridate" "forcats"
## [11] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [16] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [21] "grDevices" "utils" "datasets" "methods" "base"
##
## [[14]]
## [1] "knitr" "viridis" "viridisLite" "RColorBrewer" "effects"
## [6] "carData" "modelsummary" "fst" "infer" "lubridate"
## [11] "forcats" "stringr" "dplyr" "purrr" "readr"
## [16] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [21] "graphics" "grDevices" "utils" "datasets" "methods"
## [26] "base"
##
## [[15]]
## [1] "kableExtra" "knitr" "viridis" "viridisLite" "RColorBrewer"
## [6] "effects" "carData" "modelsummary" "fst" "infer"
## [11] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [16] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [21] "stats" "graphics" "grDevices" "utils" "datasets"
## [26] "methods" "base"
##
## [[16]]
## [1] "rmarkdown" "kableExtra" "knitr" "viridis" "viridisLite"
## [6] "RColorBrewer" "effects" "carData" "modelsummary" "fst"
## [11] "infer" "lubridate" "forcats" "stringr" "dplyr"
## [16] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [21] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [26] "datasets" "methods" "base"
##
## [[17]]
## [1] "ggridges" "rmarkdown" "kableExtra" "knitr" "viridis"
## [6] "viridisLite" "RColorBrewer" "effects" "carData" "modelsummary"
## [11] "fst" "infer" "lubridate" "forcats" "stringr"
## [16] "dplyr" "purrr" "readr" "tidyr" "tibble"
## [21] "ggplot2" "tidyverse" "stats" "graphics" "grDevices"
## [26] "utils" "datasets" "methods" "base"
##
## [[18]]
## [1] "ggridges" "rmarkdown" "kableExtra" "knitr" "viridis"
## [6] "viridisLite" "RColorBrewer" "effects" "carData" "modelsummary"
## [11] "fst" "infer" "lubridate" "forcats" "stringr"
## [16] "dplyr" "purrr" "readr" "tidyr" "tibble"
## [21] "ggplot2" "tidyverse" "stats" "graphics" "grDevices"
## [26] "utils" "datasets" "methods" "base"
##
## [[19]]
## [1] "questionr" "ggridges" "rmarkdown" "kableExtra" "knitr"
## [6] "viridis" "viridisLite" "RColorBrewer" "effects" "carData"
## [11] "modelsummary" "fst" "infer" "lubridate" "forcats"
## [16] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [21] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [26] "grDevices" "utils" "datasets" "methods" "base"ess <- read_fst("All-ESS-Data.fst")Filtering and recoding
In what follows, we will filter to France as our country of interest and re-code a number of outcomes and covariates that can be subsequently used.
Some outcomes
finand_data <- ess %>%
# Filter for France
filter(cntry == "FI") %>%
# Recoding clsprty and cleaning trstplt (will use later)
mutate(
trstplc = case_when(
trstplc == 0 ~ 1, # Keep 1 as numeric
trstplc == 2 ~ 0, # Recode 2 as numeric 0
trstplc %in% c(77, 88, 99) ~ NA_real_, # Set 7, 8, 9 as NA keeping the type numeric
TRUE ~ trstplc # Keep all other values as they are
),
trstlgl = ifelse(trstlgl %in% c(77, 88, 99), NA, trstlgl), ## Example of a trust variable
pplhlp = ifelse(pplhlp %in% c(77, 88, 99), NA, pplhlp), # Example of a satisfaction variable
medcrgvc = ifelse(medcrgvc %in% c(77, 88, 99), NA, medcrgvc), # Example of a immigration attitude variable
wrkorg = ifelse(wrkorg %in% c(7, 8, 9), NA, wrkorg), # Example of a human values variable
)#Finland_data_wrkorg
finland_data <- ess %>%
mutate(gndr_recode = case_when(
gndr == 1 ~ 'Male',
gndr == 2 ~ 'Female', # to retain other values as they are (optional)
))
# check the recoded variable table(finland_data$wrkorg_recode)
finland_data <- finland_data %>%
mutate(trstplc_recode = recode(as.character(trstplc),
'0' = "Untrustworthy",
'1' = "Moderately Trustworthy",
'2' = "Moderately Trustworthy",
'3' = "Moderately Trustworthy",
'4' = "Moderately Trustworthy",
'5' = "Moderately Trustworthy",
'6' = "Mostly Trustworthy",
'7' = "Mostly Trustworthy",
'8' = "Mostly Trustworthy",
'9' = "Mostly Trustworthy",
'10' = "Completely Trustworthy",
'77' = NA_character_,
'88' = NA_character_,
'99' = NA_character_))
# check table(finland_data$trstplc_finland)
finland_data <- finland_data %>%
mutate(pplhlp_recode = recode(as.character(pplhlp),
'0' = "Selfish",
'1' = "Moderate",
'2' = "Moderate",
'3' = "Moderate",
'4' = "Moderate",
'5' = "Moderate",
'6' = "Helpful",
'7' = "Helpful",
'8' = "Helpful",
'9' = "Helpful",
'10' = "Selfless",
'77' = NA_character_,
'88' = NA_character_,
'99' = NA_character_))
finland_data <- finland_data %>%
mutate(medcrgvc_recode = recode(as.character(medcrgvc),
'0' = "Does not apply at all",
'1' = "Moderately Free",
'2' = "Moderately Free",
'3' = "Moderately Free",
'4' = "Moderately Free",
'5' = "Moderately Free",
'6' = "Mostly Free",
'7' = "Mostly Free",
'8' = "Mostly Free",
'9' = "Mostly Free",
'10' ="COmpletely Free",
'77' = NA_character_,
'88' = NA_character_,
'99' = NA_character_))
#Check table(finland_data$medcrgvc_recode)
model1 <- lm(trstplc ~ gndr_recode, data = finland_data)
model2 <- lm(trstplc ~ gndr_recode + medcrgvc_recode, data = finland_data)
model3 <- lm(trstplc ~ gndr_recode + medcrgvc_recode + gndr_recode*medcrgvc_recode, data = finland_data)
modelsummary(
list(model1, model2, model3))| Â (1) | Â Â (2) | Â Â (3) | |
|---|---|---|---|
| (Intercept) | 7.298 | 7.589 | 7.700 |
| (0.018) | (0.057) | (0.073) | |
| gndr_recodeMale | −0.688 | −0.411 | −0.623 |
| (0.027) | (0.048) | (0.101) | |
| medcrgvc_recodeDoes not apply at all | −2.844 | −2.728 | |
| (0.133) | (0.184) | ||
| medcrgvc_recodeModerately Free | −1.484 | −1.606 | |
| (0.070) | (0.097) | ||
| medcrgvc_recodeMostly Free | −0.512 | −0.680 | |
| (0.061) | (0.087) | ||
| gndr_recodeMale × medcrgvc_recodeDoes not apply at all | −0.265 | ||
| (0.266) | |||
| gndr_recodeMale × medcrgvc_recodeModerately Free | 0.238 | ||
| (0.140) | |||
| gndr_recodeMale × medcrgvc_recodeMostly Free | 0.336 | ||
| (0.123) | |||
| Num.Obs. | 489588 | 110349 | 110349 |
| R2 | 0.001 | 0.008 | 0.008 |
| R2 Adj. | 0.001 | 0.008 | 0.008 |
| AIC | 3587570.3 | 771423.4 | 771417.9 |
| BIC | 3587603.6 | 771481.0 | 771504.4 |
| Log.Lik. | −1793782.147 | −385705.681 | −385699.955 |
| RMSE | 9.44 | 7.98 | 7.98 |
#Note: moderator vs mediators,
modelsummary(
list(model1, model2, model3),
fmt = 1,
estimate = c( "{estimate} ({std.error}){stars}",
"{estimate} ({std.error}){stars}",
"{estimate} ({std.error}){stars}"),
statistic = NULL,
coef_omit = "Intercept") | Â (1) | Â Â (2) | Â Â (3) | |
|---|---|---|---|
| gndr_recodeMale | −0.7 (0.0)*** | −0.4 (0.0)*** | −0.6 (0.1)*** |
| medcrgvc_recodeDoes not apply at all | −2.8 (0.1)*** | −2.7 (0.2)*** | |
| medcrgvc_recodeModerately Free | −1.5 (0.1)*** | −1.6 (0.1)*** | |
| medcrgvc_recodeMostly Free | −0.5 (0.1)*** | −0.7 (0.1)*** | |
| gndr_recodeMale × medcrgvc_recodeDoes not apply at all | −0.3 (0.3) | ||
| gndr_recodeMale × medcrgvc_recodeModerately Free | 0.2 (0.1)+ | ||
| gndr_recodeMale × medcrgvc_recodeMostly Free | 0.3 (0.1)** | ||
| Num.Obs. | 489588 | 110349 | 110349 |
| R2 | 0.001 | 0.008 | 0.008 |
| R2 Adj. | 0.001 | 0.008 | 0.008 |
| AIC | 3587570.3 | 771423.4 | 771417.9 |
| BIC | 3587603.6 | 771481.0 | 771504.4 |
| Log.Lik. | −1793782.147 | −385705.681 | −385699.955 |
| RMSE | 9.44 | 7.98 | 7.98 |
### Seperate here
# We use the "effect" function to call the interaction term from model 3. We assign the results to "effect_plot"
interaction_plot <- effect("gndr_recode*medcrgvc_recode", model3, na.rm=TRUE)
# We plot here the interaction term
plot(interaction_plot,
main="Interaction effect Female and Male Gender x level of media criticism",
xlab="Male and Female Gender",
ylab="media criticism")
### Plotting
plot <- expand.grid(
pplhlp_recode = unique(finland_data$pplhlp_recode), # Using cohort here
medcrgvc_recode = unique(finland_data$medcrgvc_recode)
)
model4 <- lm(trstplc ~ pplhlp_recode + medcrgvc_recode + pplhlp_recode*medcrgvc_recode, data = finland_data)
plot$predicted <- predict(model4, newdata = plot)
# Aggregate data by cohort and educ_level
avg_data <- finland_data %>%
group_by(pplhlp_recode, medcrgvc_recode) %>%
summarize(avg_trstplc = mean(trstplc, na.rm = TRUE))## `summarise()` has grouped output by 'pplhlp_recode'. You can override using the
## `.groups` argument.# Plotting
ggplot(avg_data, aes(x = pplhlp_recode, y = avg_trstplc, color = medcrgvc_recode)) +
geom_point(size = 0.75) +
geom_line(data = plot, aes(y = predicted, group = medcrgvc_recode)) +
scale_y_continuous(limits = c(3.5, 7), breaks = scales::breaks_width(0.5)) +
scale_color_discrete(breaks = c("Yes", "No")) +
theme(legend.position = "bottom") +
labs(
x = "Level of helpfulness",
y = "trust in police",
color = "level of media criticism"
)## Warning: Removed 12 rows containing missing values (`geom_point()`).## Warning: Removed 14 rows containing missing values (`geom_line()`).
```