Lei_Amy_Homework_4
2024-03-18
Setting up environment
# List of packages
packages <- c("tidyverse", "infer", "fst", "modelsummary", "effects", "survey", "interactions","tidymodels") # 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.4 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.4.4 ✔ 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: grid
##
## Loading required package: Matrix
##
##
## Attaching package: 'Matrix'
##
##
## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack
##
##
## Loading required package: survival
##
##
## Attaching package: 'survey'
##
##
## The following object is masked from 'package:graphics':
##
## dotchart
##
##
## Registered S3 methods overwritten by 'broom':
## method from
## tidy.glht jtools
## tidy.summary.glht jtools
##
## ── Attaching packages ────────────────────────────────────── tidymodels 1.1.1 ──
##
## ✔ broom 1.0.5 ✔ rsample 1.2.0
## ✔ dials 1.2.1 ✔ tune 1.1.2
## ✔ modeldata 1.3.0 ✔ workflows 1.1.4
## ✔ parsnip 1.2.0 ✔ workflowsets 1.0.1
## ✔ recipes 1.0.10 ✔ yardstick 1.3.0
##
## ── Conflicts ───────────────────────────────────────── tidymodels_conflicts() ──
## ✖ scales::discard() masks purrr::discard()
## ✖ Matrix::expand() masks tidyr::expand()
## ✖ dplyr::filter() masks stats::filter()
## ✖ recipes::fixed() masks stringr::fixed()
## ✖ dplyr::lag() masks stats::lag()
## ✖ Matrix::pack() masks tidyr::pack()
## ✖ yardstick::spec() masks readr::spec()
## ✖ recipes::step() masks stats::step()
## ✖ Matrix::unpack() masks tidyr::unpack()
## ✖ recipes::update() masks Matrix::update(), stats::update()
## • Search for functions across packages at https://www.tidymodels.org/find/## [[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] "survey" "survival" "Matrix" "grid" "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"
##
## [[7]]
## [1] "interactions" "survey" "survival" "Matrix" "grid"
## [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"
##
## [[8]]
## [1] "yardstick" "workflowsets" "workflows" "tune" "rsample"
## [6] "recipes" "parsnip" "modeldata" "dials" "scales"
## [11] "broom" "tidymodels" "interactions" "survey" "survival"
## [16] "Matrix" "grid" "effects" "carData" "modelsummary"
## [21] "fst" "infer" "lubridate" "forcats" "stringr"
## [26] "dplyr" "purrr" "readr" "tidyr" "tibble"
## [31] "ggplot2" "tidyverse" "stats" "graphics" "grDevices"
## [36] "utils" "datasets" "methods" "base"Task 1
Graph the histogram distribution using the populist scale (as we did in the tutorial) for the country of Sweden. What do you note?
sweden_data <- read.fst("sweden_data.fst")sweden_data$trstplt <- ifelse(sweden_data$trstplt > 10, NA, sweden_data$trstplt)
sweden_data$trstprl <- ifelse(sweden_data$trstprl > 10, NA, sweden_data$trstprl)
sweden_data$trstprt <- ifelse(sweden_data$trstprt > 10, NA, sweden_data$trstprt)
sweden_data$trust <- scales::rescale(sweden_data$trstplt + sweden_data$trstprl + sweden_data$trstprt, na.rm = TRUE, to = c(0, 100))
sweden_data$populist <- scales::rescale(sweden_data$trust, na.rm = TRUE, to=c(100,0))ggplot(sweden_data, aes(x = populist)) +
geom_histogram(bins = 30, fill = "blue", color = "black") +
theme_minimal() +
labs(title = "Distribution of Populist Scale for Sweden",
x = "Populist Scale",
y = "Count")## Warning: Removed 2503 rows containing non-finite values (`stat_bin()`).
The distribution depicted by the histogram indicates a minor skew to the
right, with the most frequent value (mode) being below 50. This suggests
a low degree of trust among the majority of the Swedish population
towards their government, including parliaments, politicians, and
political parties. Additionally, the data shows that the segments
representing 0 to 25 and 75 to 100 on the populist scale have lesser
instances compared to the 25 to 75 range. Specifically, the 25 to 50
range on the scale records a higher frequency than the 50 to 75
segment.
Task 2
Run a linear regression with populist attitudes as the outcome, educational attainment (recoded as a binary “BA or more” and “No BA”), and using the Swedish data. Print the intercept and coefficients only and interpret. Then do a tidy model table displaying the p-value (with digits = 3) and interpret.
sweden_data <- read.fst("sweden_data.fst")# Populist scale (Norris/Inglehart use low levels of trust in parliaments, politicians and parties as an indicator of populism)
# Note: Schafer does it 2 ways, we will focus now on the most straightforward
# Setting values greater than 10 to NA
sweden_data$trstplt <- ifelse(sweden_data$trstplt > 10, NA, sweden_data$trstplt)
sweden_data$trstprl <- ifelse(sweden_data$trstprl > 10, NA, sweden_data$trstprl)
sweden_data$trstprt <- ifelse(sweden_data$trstprt > 10, NA, sweden_data$trstprt)
# Creating and rescaling the trust variable
sweden_data$trust <- scales::rescale(sweden_data$trstplt + sweden_data$trstprl + sweden_data$trstprt, na.rm = TRUE, to = c(0, 100))
# Here's how you would "flip it" because populist is conceived by N + I as 'anti-trust'
# As Schafer explains there are some issues with that conceptualization, but N + I is also widely applied
sweden_data$populist <- scales::rescale(sweden_data$trust, na.rm = TRUE, to=c(100,0))sweden_data <- sweden_data %>%
mutate(
# Recoding education
educ.ba = case_when(
essround < 5 & edulvla == 5 ~ "BA or more",
essround >= 5 & edulvlb > 600 ~ "BA or more",
TRUE ~ "No BA"
),
# Handle NAs for education levels
edulvla = ifelse(edulvla %in% c(77, 88, 99), NA_integer_, edulvla),
edulvlb = ifelse(edulvlb %in% c(5555, 7777, 8888), NA_integer_, edulvlb),
# Explicitly making 'No BA' the reference category
educ.ba = factor(educ.ba, levels = c("No BA", "BA or more")),
)
sweden_pop_model <- lm(populist ~ educ.ba, data = sweden_data)
coefficients_sweden <- coef(sweden_pop_model)
print(coefficients_sweden)## (Intercept) educ.baBA or more
## 51.358680 -7.868121tidy(sweden_pop_model)## # A tibble: 2 × 5
## term estimate std.error statistic p.value
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 51.4 0.181 283. 0
## 2 educ.baBA or more -7.87 0.350 -22.5 2.58e-110tidy(sweden_pop_model) |>
knitr::kable(digits = 3)| term | estimate | std.error | statistic | p.value |
|---|---|---|---|---|
| (Intercept) | 51.359 | 0.181 | 283.363 | 0 |
| educ.baBA or more | -7.868 | 0.350 | -22.494 | 0 |
summary(sweden_pop_model)##
## Call:
## lm(formula = populist ~ educ.ba, data = sweden_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -51.359 -14.692 -1.359 11.975 56.509
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 51.3587 0.1812 283.36 <2e-16 ***
## educ.baBA or more -7.8681 0.3498 -22.49 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.43 on 15711 degrees of freedom
## (2503 observations deleted due to missingness)
## Multiple R-squared: 0.0312, Adjusted R-squared: 0.03114
## F-statistic: 506 on 1 and 15711 DF, p-value: < 2.2e-16The intercept value of 51.4 in this model provides the baseline average score on the populist scale for the reference demographic in Sweden—individuals without a Bachelor’s degree (“No BA”). This suggests that those with lower educational attainment have a propensity towards populist attitudes, with the average score reflecting moderate populist sentiment. Higher scores on this scale are indicative of more pronounced populist views.
The coefficient of -7.9 for individuals possessing a Bachelor’s degree or higher reflects a clear divergence from this baseline. Specifically, those with higher educational attainment in Sweden score, on average, 7.9 points lower on the populist scale compared to their counterparts without a Bachelor’s degree. This finding suggests that higher education is associated with less populist leanings.
The p-value, registered at 0 for both educational groups, underscores the statistical significance of this relationship. The absence of a p-value above the alpha threshold implies a meaningful correlation between levels of education and populist attitudes within the Swedish context. The results point to a discernible pattern where educational attainment influences populist sentiments, with higher education seemingly contributing to a decrease in such attitudes.
Task 3
Now using the Italian dataset and the authoritarian values scale (as we did in the tutorial), graph the average by survey year. Interpret.
df <- read.fst("italy_data.fst")
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 ))
# Now you can calculate 'schwartzauth' after the NA recoding
df$auth <- scales::rescale(df$behave +
df$secure +
df$safety +
df$tradition +
df$rules, to=c(0,100), na.rm=TRUE)
df <- df %>% filter(!is.na(auth))
table(df$secure)##
## 1 2 3 4 5 6
## 54 139 626 1856 2876 2881table(df$auth)##
## 0 12 16 20 24 28 32 36 40 44 48 52 56 60 64 68 72 76 80 84
## 7 4 4 9 13 20 27 57 113 155 163 323 367 596 647 789 823 897 992 733
## 88 92 96 100
## 632 465 299 297df$year <- NA
replacements <- c(2002, 2004, 2006, 2008, 2010, 2012, 2014, 2016, 2018, 2020)
for(i in 1:10){
df$year[df$essround == i] <- replacements[i]
}
auth_avg <- df %>%
group_by(year) %>%
summarize(auth_avg = mean(auth, na.rm = TRUE))
plot_auth <- ggplot(auth_avg, aes(x = year, y = auth_avg)) +
geom_point(aes(color = auth_avg), alpha = 1.0) +
geom_line(aes(group = 1), color = "blue", linetype = "dashed") +
labs(title = "Authoritarian Scale Average by Year (Italy)",
x = "Survey Year",
y = "Authoritarian Scale Average") +
theme_minimal() +
theme(legend.position = "none") +
scale_y_continuous(limits = c(0, 100))
print(plot_auth)
The trend of the authoritarian scale average in Italy over the years presents a stable pattern from 2012 through 2020, indicating minimal fluctuations. Initially, there’s a minor dip in the scale from 75 to approximately 73 over the four years leading up to 2016. However, this downward trend reverses in the subsequent period, from 2016 to 2020, where the scale gradually ascends back to a level close to 75. This consistency in the authoritarian scale suggests that the average level remains relatively high, hovering around the 75 mark throughout the observed period.
##Task 4
Now run a linear regression model, using the modelsummary package, with authoritarian values as the outcome, and generations as the predictor/potential explanatory variable, and using the Italian data again. Rename the coefficients using the coef_rename function. Interpret the coefficients, whether they are statistically significant, and the adjusted R-squared.
model_data <- df %>%
mutate(
# Recoding cohort variable, setting years before 1930 and after 2000 to NA
cohort = ifelse(yrbrn < 1930 | yrbrn > 2000, NA, yrbrn),
# Recoding generational cohorts based on year of birth
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(cohort) # Keeping other values as character if they do not fit the ranges
),
# Converting cohort to a factor with labels
gen = factor(gen,
levels = c("1", "2", "3", "4"),
labels = c("Interwar", "Baby Boomers", "Gen X", "Millennials"))
)
table(model_data$gen)##
## Interwar Baby Boomers Gen X Millennials
## 1026 2580 2214 1811model <- lm(auth ~ gen, data = model_data)
modelsummary(model,
fmt = 1,
estimate = c("{estimate} ({std.error}){stars}"),
statistic = NULL,
coef_rename = c("genBaby Boomers" = "Baby Boomers", "genGen X" = "Gen X", "genMillennials" = "Millennials"),
title = 'Table1. Regression model for authoritarian attitudes')| (1) | |
|---|---|
| (Intercept) | 75.3 (0.5)*** |
| Baby Boomers | −1.1 (0.6)* |
| Gen X | −2.1 (0.6)*** |
| Millennials | −3.6 (0.6)*** |
| Num.Obs. | 7631 |
| R2 | 0.006 |
| R2 Adj. | 0.006 |
| AIC | 62934.0 |
| BIC | 62968.7 |
| Log.Lik. | −31462.009 |
| RMSE | 14.94 |
The baseline authoritarian scale score, set at 75.3, serves as the predicted average score in the absence of any influence from the examined predictors. This figure shifts when accounting for generational differences, with Baby Boomers displaying authoritarian scale scores that are, on average, 1.1 points below that of the reference category, assuming other factors remain unchanged. For individuals belonging to Generation X, this difference widens to -2.1 points, indicating a further decrease in authoritarian scale scores relative to the reference group with all other variables held constant. Millennials exhibit an even larger deviation, with their scores averaging 3.6 points below the reference group’s score, suggesting a significant generational shift in authoritarian scale scores when controlling for other influences.
##Task 5
Now visualize, using the modelplot() function from the modelsummary package, coefficient estimates and 95% confidence intervals for the following three models:
Model 1 = populist attitudes scale as the outcome; left/right as the single predictor (omitting moderates as we did in the tutorial), and using data for Greece.
Model 2 = populist attitudes scale as the outcome; gender as the single predictor, and using data for Greece.
Model 3 = populist attitudes scale as the outcome; generations as the single predictor (using the four generational category recode we did in the tutorial), and using data for Greece. Interpret.
greece_data <- read.fst("greece_data.fst")greece_data$trstplt <- ifelse(greece_data$trstplt > 10, NA, greece_data$trstplt)
greece_data$trstprl <- ifelse(greece_data$trstprl > 10, NA, greece_data$trstprl)
greece_data$trstprt <- ifelse(greece_data$trstprt > 10, NA, greece_data$trstprt)
greece_data$trust <- scales::rescale(greece_data$trstplt + greece_data$trstprl + greece_data$trstprt, na.rm = TRUE, to = c(0, 100))
greece_data$populist <- scales::rescale(greece_data$trust, na.rm = TRUE, to=c(100,0))model_data <- greece_data %>%
mutate(
# Recoding gender
gndr = case_when(
gndr == 1 ~ "Male",
gndr == 2 ~ "Female",
gndr == 9 ~ NA_character_,
TRUE ~ as.character(gndr)
),
lrscale = case_when(
lrscale %in% 0:4 ~ "Left", # This time: Values 0 to 3 in 'lrscale' are categorized as "Left"
lrscale %in% 6:10 ~ "Right", # This time: Values 6 to 10 in 'lrscale' are categorized as "Right"
lrscale %in% c(77, 88, 99) ~ NA_character_, # Values 77, 88, 99 in 'lrscale' are set as NA
TRUE ~ NA_character_ # Any other values not explicitly matched above
),
# Recoding cohort variable, setting years before 1930 and after 2000 to NA
cohort = ifelse(yrbrn < 1930 | yrbrn > 2000, NA, yrbrn),
# Recoding generational cohorts based on year of birth
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(cohort) # Keeping other values as character if they do not fit the ranges
),
# Converting cohort to a factor with labels
gen = factor(gen,
levels = c("1", "2", "3", "4"),
labels = c("Interwar", "Baby Boomers", "Gen X", "Millennials"))
)
table(model_data$gen)##
## Interwar Baby Boomers Gen X Millennials
## 2978 3484 3498 2330odel1 <- lm(populist ~ lrscale, data = model_data)
model2 <- lm(populist ~ gndr, data = model_data)
model3 <- lm(populist ~ gen, data = model_data)models <- list(
"Model 1" = lm(populist ~ lrscale, data = model_data),
"Model 2" = lm(populist ~ gndr, data = model_data),
"Model 3" = lm(populist ~ gen, data = model_data))
modelsummary(models, fmt = 1,
estimate = c("{estimate} ({std.error}){stars}",
"{estimate} ({std.error}){stars}",
"{estimate} ({std.error}){stars}"),
statistic = NULL)| Model 1 | Model 2 | Model 3 | |
|---|---|---|---|
| (Intercept) | 72.0 (0.5)*** | 70.0 (0.3)*** | 65.8 (0.5)*** |
| lrscaleRight | −8.8 (0.6)*** | ||
| gndrMale | −1.0 (0.5)* | ||
| genBaby Boomers | 3.5 (0.7)*** | ||
| genGen X | 5.6 (0.7)*** | ||
| genMillennials | 5.8 (0.7)*** | ||
| Num.Obs. | 4936 | 9823 | 9562 |
| R2 | 0.037 | 0.000 | 0.009 |
| R2 Adj. | 0.036 | 0.000 | 0.009 |
| AIC | 44723.8 | 88873.2 | 86473.9 |
| BIC | 44743.3 | 88894.8 | 86509.7 |
| Log.Lik. | −22358.894 | −44433.590 | −43231.944 |
| RMSE | 22.44 | 22.30 | 22.25 |
modelplot(models, coef_omit = 'Interc')
The graphical representation illustrates a stratified analysis of generational impact on authoritarian values, set against a reference category. The data points for the Baby Boomer generation indicate a notable uptick of 3.5 points on the authoritarian values scale, signifying a statistical departure from the values held by the Interwar generation. This divergence is marked by substantial statistical significance, as emphasized by the annotation of three asterisks, reflecting a p-value below the conventional threshold of 0.05.
Firstly, for Generation X, we observe an accentuated departure from the Interwar baseline, with an average increase of 5.3 points in authoritarian values. The presence of three asterisks here again denotes a level of statistical significance that falls well below the 0.05 alpha threshold, underscoring the robustness of the observed generational shift.
Millennials exhibit the most substantial rise in authoritarian values, standing at 5.8 points above the reference group. The statistical significance of this increase is considerable, as indicated by the same notation of three asterisks, suggesting a strong deviation in authoritarian values among Millennials compared to the Interwar generation.
The model’s adjusted R-squared value of 0.009 suggests that while the generational category does have a statistically significant effect on the authoritarian values scale in Greece, it accounts for just under 1% of the variation. This modest explanatory power indicates that additional factors, beyond generational categorization, likely contribute to shaping authoritarian values within the population.
In terms of gender, the coefficient associated with being male is not statistically significant, suggesting that within this model, gender does not appear to have a discernible effect on the authoritarian values scale when compared to the female category. Conversely, the political orientation model indicates a substantial statistical significance; those aligned with the right, as opposed to the left (reference category), have scores that are, on average, 8.8 points lower on the authoritarian scale. Despite this significant finding, the model accounts for 3.6% of the variance in authoritarian values, pointing again to the existence of other influential factors not captured in this analysis. This model’s relative explanatory strength is greater than that of the generational model, yet it still leaves the majority of variance unaccounted for, suggesting a complex interplay of variables in determining authoritarian values.