SOC222 REPORT FINDINGS
Fatima Nawaz
2024-03-29
Reearch Q: “How does political ideology relate to trust in science in Iceland?”
Variables
IV (explanatory): Political Ideology (lrscale)
DV (outcome): Trust in Scientists (trstsci)
Interaction Variable: Education (educ.ba)
Loading Packages
# List of packages
packages <- c("tidyverse", "infer", "fst", "modelsummary", "effects", "survey", "performance", "flextable", "broom", "effects", "scales", "ggridges", "viridis", "knitr", "kableExtra", "rmarkdown", "questionr")
# 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.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.0 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ✔ 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
## Version 2.0.0 of `modelsummary`, to be released soon, will introduce a
## breaking change: The default table-drawing package will be `tinytable`
## instead of `kableExtra`. All currently supported table-drawing packages
## will continue to be supported for the foreseeable future, including
## `kableExtra`, `gt`, `huxtable`, `flextable, and `DT`.
##
## You can always call the `config_modelsummary()` function to change the
## default table-drawing package in persistent fashion. To try `tinytable`
## now:
##
## config_modelsummary(factory_default = 'tinytable')
##
## To set the default back to `kableExtra`:
##
## config_modelsummary(factory_default = 'kableExtra')
##
## 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
##
##
##
## Attaching package: 'flextable'
##
##
## The following object is masked from 'package:purrr':
##
## compose
##
##
##
## Attaching package: 'scales'
##
##
## The following object is masked from 'package:purrr':
##
## discard
##
##
## The following object is masked from 'package:readr':
##
## col_factor
##
##
## Loading required package: viridisLite
##
##
## Attaching package: 'viridis'
##
##
## The following object is masked from 'package:scales':
##
## viridis_pal
##
##
##
## Attaching package: 'kableExtra'
##
##
## The following objects are masked from 'package:flextable':
##
## as_image, footnote
##
##
## 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] "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] "performance" "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] "flextable" "performance" "survey" "survival" "Matrix"
## [6] "grid" "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"
##
## [[9]]
## [1] "broom" "flextable" "performance" "survey" "survival"
## [6] "Matrix" "grid" "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"
##
## [[10]]
## [1] "broom" "flextable" "performance" "survey" "survival"
## [6] "Matrix" "grid" "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"
##
## [[11]]
## [1] "scales" "broom" "flextable" "performance" "survey"
## [6] "survival" "Matrix" "grid" "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"
##
## [[12]]
## [1] "ggridges" "scales" "broom" "flextable" "performance"
## [6] "survey" "survival" "Matrix" "grid" "effects"
## [11] "carData" "modelsummary" "fst" "infer" "lubridate"
## [16] "forcats" "stringr" "dplyr" "purrr" "readr"
## [21] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [26] "graphics" "grDevices" "utils" "datasets" "methods"
## [31] "base"
##
## [[13]]
## [1] "viridis" "viridisLite" "ggridges" "scales" "broom"
## [6] "flextable" "performance" "survey" "survival" "Matrix"
## [11] "grid" "effects" "carData" "modelsummary" "fst"
## [16] "infer" "lubridate" "forcats" "stringr" "dplyr"
## [21] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [26] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [31] "datasets" "methods" "base"
##
## [[14]]
## [1] "knitr" "viridis" "viridisLite" "ggridges" "scales"
## [6] "broom" "flextable" "performance" "survey" "survival"
## [11] "Matrix" "grid" "effects" "carData" "modelsummary"
## [16] "fst" "infer" "lubridate" "forcats" "stringr"
## [21] "dplyr" "purrr" "readr" "tidyr" "tibble"
## [26] "ggplot2" "tidyverse" "stats" "graphics" "grDevices"
## [31] "utils" "datasets" "methods" "base"
##
## [[15]]
## [1] "kableExtra" "knitr" "viridis" "viridisLite" "ggridges"
## [6] "scales" "broom" "flextable" "performance" "survey"
## [11] "survival" "Matrix" "grid" "effects" "carData"
## [16] "modelsummary" "fst" "infer" "lubridate" "forcats"
## [21] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [26] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [31] "grDevices" "utils" "datasets" "methods" "base"
##
## [[16]]
## [1] "rmarkdown" "kableExtra" "knitr" "viridis" "viridisLite"
## [6] "ggridges" "scales" "broom" "flextable" "performance"
## [11] "survey" "survival" "Matrix" "grid" "effects"
## [16] "carData" "modelsummary" "fst" "infer" "lubridate"
## [21] "forcats" "stringr" "dplyr" "purrr" "readr"
## [26] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [31] "graphics" "grDevices" "utils" "datasets" "methods"
## [36] "base"
##
## [[17]]
## [1] "questionr" "rmarkdown" "kableExtra" "knitr" "viridis"
## [6] "viridisLite" "ggridges" "scales" "broom" "flextable"
## [11] "performance" "survey" "survival" "Matrix" "grid"
## [16] "effects" "carData" "modelsummary" "fst" "infer"
## [21] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [26] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [31] "stats" "graphics" "grDevices" "utils" "datasets"
## [36] "methods" "base"Loading and filtering data for ICELAND = our country of interest
ess <- read_fst("All-ESS-Data.fst")
iceland_data <- ess |>
filter(cntry == "IS") |> # "IS" = Iceland
select(edulvla, edulvlb, trstsci, lrscale, essround, pweight, dweight) # These are all the variables we selected out of the entire ESS datasetCoding ‘lrscale’ as our explanatory variable and ‘educ.ba’ as our interaction variable (additional predictor)
iceland_data_clean <- iceland_data |>
mutate(
lrscale = case_when(
lrscale == 0 ~ "Left",
lrscale == 1 ~ "Left",
lrscale == 2 ~ "Left",
lrscale == 3 ~ "Left",
lrscale == 4 ~ "Centrist",
lrscale == 5 ~ "Centrist",
lrscale == 6 ~ "Centrist",
lrscale == 7 ~ "Right",
lrscale == 8 ~ "Right",
lrscale == 9 ~ "Right",
lrscale == 10 ~ "Right",
lrscale %in% c(77, 88, 99) ~ NA_character_,
TRUE ~ as.character(lrscale)
),
educ.ba = case_when(
essround < 5 & edulvla == 5 ~ "BA or more",
essround >= 5 & edulvlb > 600 ~ "BA or more",
TRUE ~ "No BA"
)
) |>
mutate(
edulvla = ifelse(edulvla %in% c(77,88), NA_integer_, edulvla),
edulvlb = ifelse(edulvlb %in% c(5555, 7777, 8888), NA_integer_, edulvlb),
educ.ba = as.character(educ.ba)
)
table(iceland_data_clean$educ.ba)##
## BA or more No BA
## 1351 2624table(iceland_data_clean$lrscale)##
## Centrist Left Right
## 1796 899 1033Re-naming variables
iceland_data_clean <- iceland_data_clean |>
mutate(Education = educ.ba,
Trust_in_Scientists = trstsci,
Ideology = lrscale) |>
mutate(across(c("Trust_in_Scientists"),
~ na_if(.x, 99) |> na_if(88) |> na_if(77))) Double-checking if ‘NAs’ are clean
unique(iceland_data_clean$Education)## [1] "No BA" "BA or more"unique(iceland_data_clean$Ideology)## [1] NA "Centrist" "Left" "Right"unique(iceland_data_clean$Trust_in_Scientists)## [1] NA 8 10 7 5 6 9 0 3 4 2iceland_data_clean <- iceland_data_clean %>% filter(!is.na(Ideology))
iceland_data_clean <- iceland_data_clean %>% filter(!is.na(Trust_in_Scientists))Calculating Weights
iceland_data_clean$weight <- iceland_data_clean$dweight * iceland_data_clean$pweight # creating weight column based on multiplication of d and p weight
survey_design <- svydesign(ids = ~1, data = iceland_data_clean, weights = ~weight) # telling R how to interpret dataset - identifying different parts of the data, R knows what a weight means, so just have to tell them what's whatNEXT: INFER PACKAGE STATISTICAL ASSOCIATION TEST
To assess whether there is a significant association between political ideology and trust in scientists that goes beyond random chance.
We are doing ANOVA because our explanatory variable is three categories and our response variable is numeric, so according to the table in our lecture that corresponds to the ANOVA test statistic
Step 1: Calculate the test statistic of your sample
test_stat_anova <- iceland_data_clean |>
specify(explanatory = Ideology, # change variable name for explanatory variable
response = Trust_in_Scientists) |> # change variable name for outcome of interest
hypothesize(null = "independence") |> # store data until we do a test
calculate(stat = "F") # replace in between quotation marks appropriate test statistic
print(test_stat_anova$stat)## [1] 9.920009Interpreting Value
9.92 is the observed F-test statistic that we will compare to the null distribution to determine whether there is significance between political ideology and trust in scientists in Iceland.
Step 2: Simulate the null distribution
null_distribution_anova <- iceland_data_clean |>
specify(explanatory = Ideology,
response = Trust_in_Scientists) |>
hypothesize(null = "independence") |>
generate(reps = 1000, type = "permute") |>
calculate(stat = "F") # 1000 different times it calculates ANOVA statisticInterpreting Null Simulation
This is creating a simulation of the distribution of expected values of ‘Ideology’ and ’Trust_in_Scientists” when they are cross-tabulated, under the hypothetical assumption that they are independent. This means assuming there is no difference in political ideology and trust in scientists for Iceland.
Step 3: Calculate the p-value of your sample
p_val_anova <- null_distribution_anova |>
get_pvalue(obs_stat = test_stat_anova, direction = "greater") ## Warning: Please be cautious in reporting a p-value of 0. This result is an approximation
## based on the number of `reps` chosen in the `generate()` step.
## ℹ See `get_p_value()` (`?infer::get_p_value()`) for more information.# given that 13 = calculated observed value how many times does the expected test statistics exceed the critical value
p_val_anova## # A tibble: 1 × 1
## p_value
## <dbl>
## 1 0Interpreting P-Value
The p-value is very close to 0, which indicates an almost 0% chance that the observed F-test statistic of 9.92 will be captured in the null distribution, suggesting that there is a significant association between political ideology and trust in scientists in Iceland.
Step 4: Visualize without confidence intervals
null_distribution_anova |>
visualize() +
shade_p_value(obs_stat = test_stat_anova, direction = "greater") # add on observed statistic and direction - two sided test## Warning in (function (mapping = NULL, data = NULL, stat = "identity", position = "identity", : All aesthetics have length 1, but the data has 1000 rows.
## ℹ Did you mean to use `annotate()`?## Warning in min(diff(unique_loc)): no non-missing arguments to min; returning
## Inf
null_distribution_anova## Response: Trust_in_Scientists (numeric)
## Explanatory: Ideology (factor)
## Null Hypothesis: independence
## # A tibble: 1,000 × 2
## replicate stat
## <int> <dbl>
## 1 1 0.384
## 2 2 0.723
## 3 3 1.10
## 4 4 1.69
## 5 5 0.697
## 6 6 0.617
## 7 7 0.597
## 8 8 3.36
## 9 9 1.50
## 10 10 1.78
## # ℹ 990 more rowsInterpretation
According to the null distribution, the observed test statistic of 9.92 (the red line) falls outside of the null distribution. This suggests there is significant evidence to reject the null hypothesis under assumptions of frequentist hypothesis testing because if the null hypothesis is true, our observed test result is significantly different from that expectation.
Thus, it suggests that trust in scientists in Iceland is can be explained by political ideology.
Step 4: Visualize with confidence intervals
conf_int_anova <- null_distribution_anova |>
get_confidence_interval(level = 0.95, type = "percentile") # showing confidence interval - if red line was in green boundary - cannot reject null
null_distribution_anova |>
visualize() +
shade_p_value(obs_stat = test_stat_anova, direction = "greater") +
shade_confidence_interval(endpoints = conf_int_anova)## Warning in (function (mapping = NULL, data = NULL, stat = "identity", position = "identity", : All aesthetics have length 1, but the data has 1000 rows.
## ℹ Did you mean to use `annotate()`?## Warning in min(diff(unique_loc)): no non-missing arguments to min; returning
## Inf
null_distribution_anova## Response: Trust_in_Scientists (numeric)
## Explanatory: Ideology (factor)
## Null Hypothesis: independence
## # A tibble: 1,000 × 2
## replicate stat
## <int> <dbl>
## 1 1 0.384
## 2 2 0.723
## 3 3 1.10
## 4 4 1.69
## 5 5 0.697
## 6 6 0.617
## 7 7 0.597
## 8 8 3.36
## 9 9 1.50
## 10 10 1.78
## # ℹ 990 more rowsInterpreting Output
The confidence interval showcases the range of not rejecting the null in the simulated null distribution. Thus, anything outside of the two boundaries of the confidence interval indicates the rejection region of the ANOVA test. Since our observed F-test statistic of 9.92 is outside of the confidence interval range, it means it is within the rejection region, so the null hypothesis is rejected.
NEXT: EXPORTING DESCRIPTIVE TABLES
Categorical Descriptive Table
table1 <- datasummary_skim(iceland_data_clean |> dplyr::select(Education, Trust_in_Scientists, Ideology), type = "categorical",
title = "Table 1. Categorical Descriptive Statistics for Main Variables", output = "flextable")table1
|
| N | % |
|---|---|---|---|
Education | BA or more | 328 | 39.0 |
No BA | 514 | 61.0 | |
Ideology | Centrist | 376 | 44.7 |
Left | 223 | 26.5 | |
Right | 243 | 28.9 |
set_flextable_defaults(fonts_ignore=TRUE)
flextable::save_as_docx(table1, path = "Table1.docx",
width = 7.0, height = 7.0)Numeric Descriptive Table
table1b <- datasummary_skim(iceland_data_clean |> dplyr::select(Education, Trust_in_Scientists, Ideology), title = "Table 1. Numeric Descriptive Statistics for Main Variables", output = "flextable")## Warning: The histogram argument is only supported for (a) output types "default",
## "html", "kableExtra", or "gt"; (b) writing to file paths with extensions
## ".html", ".jpg", or ".png"; and (c) Rmarkdown, knitr or Quarto documents
## compiled to PDF (via kableExtra) or HTML (via kableExtra or gt). Use
## `histogram=FALSE` to silence this warning.table1b
| Unique | Missing Pct. | Mean | SD | Min | Median | Max |
|---|---|---|---|---|---|---|---|
Trust_in_Scientists | 10 | 0 | 8.1 | 1.5 | 0.0 | 8.0 | 10.0 |
set_flextable_defaults(fonts_ignore=TRUE)
flextable::save_as_docx(table1b, path = "Table1b.docx",
width = 7.0, height = 7.0)NEXT: MULTIPLE LINEAR REGRESSION MODELS
Let’s run our linear regression models to assess the magnitude and direction of our signfiicant association between political ideology and trust in scientists
Let’s run our models
Creating Models
model1 <- lm(Trust_in_Scientists ~ Ideology, data = iceland_data_clean, weights = weight) # outcome + predictor #1
model2 <- lm(Trust_in_Scientists ~ Ideology + Education, data = iceland_data_clean, weights = weight) # outcome + predictor #1 + predictor #2
model3 <- lm(Trust_in_Scientists ~ Ideology + Education + Ideology*Education, data = iceland_data_clean, weights = weight) # outcome + predictor #1 + predictor #2 and interactionProducing modelsummary table
models <- list(
"Model 1" = lm(Trust_in_Scientists ~ Ideology, data = iceland_data_clean, weights = weight),
"Model 2" = lm(Trust_in_Scientists ~ Ideology + Education, data = iceland_data_clean, weights = weight),
"Model 3" = lm(Trust_in_Scientists ~ Ideology + Education + Ideology*Education, data = iceland_data_clean, weights = weight))
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) | 7.9 (0.1)*** | 8.1 (0.1)*** | 8.3 (0.1)*** |
| IdeologyLeft | 0.5 (0.1)*** | 0.5 (0.1)*** | 0.3 (0.2) |
| IdeologyRight | 0.1 (0.1) | 0.1 (0.1) | −0.2 (0.2) |
| EducationNo BA | −0.4 (0.1)*** | −0.6 (0.2)*** | |
| IdeologyLeft × EducationNo BA | 0.3 (0.3) | ||
| IdeologyRight × EducationNo BA | 0.4 (0.3)+ | ||
| Num.Obs. | 842 | 842 | 842 |
| R2 | 0.023 | 0.036 | 0.040 |
| R2 Adj. | 0.020 | 0.033 | 0.034 |
| AIC | 3046.5 | 3036.7 | 3037.4 |
| BIC | 3065.4 | 3060.4 | 3070.6 |
| Log.Lik. | −1519.244 | −1513.374 | −1511.715 |
| RMSE | 1.47 | 1.46 | 1.45 |
Visualizing coefficients
modelplot(models, coef_omit = 'Interc')
modelsummary(models, fmt = 1,
estimate = c("{estimate} ({std.error}){stars}",
"{estimate} ({std.error}){stars}",
"{estimate} ({std.error}){stars}"),
statistic = NULL,
coef_rename = c("IdeologyLeft" = "Left", "IdeologyRight" = "Right", "EducationNo BA" = "No BA"),
title = 'Table 2. Regression models for Trust In Scientists')| Model 1 | Model 2 | Model 3 | |
|---|---|---|---|
| (Intercept) | 7.9 (0.1)*** | 8.1 (0.1)*** | 8.3 (0.1)*** |
| Left | 0.5 (0.1)*** | 0.5 (0.1)*** | 0.3 (0.2) |
| Right | 0.1 (0.1) | 0.1 (0.1) | −0.2 (0.2) |
| No BA | −0.4 (0.1)*** | −0.6 (0.2)*** | |
| Left:No BA | 0.3 (0.3) | ||
| Right:No BA | 0.4 (0.3)+ | ||
| Num.Obs. | 842 | 842 | 842 |
| R2 | 0.023 | 0.036 | 0.040 |
| R2 Adj. | 0.020 | 0.033 | 0.034 |
| AIC | 3046.5 | 3036.7 | 3037.4 |
| BIC | 3065.4 | 3060.4 | 3070.6 |
| Log.Lik. | −1519.244 | −1513.374 | −1511.715 |
| RMSE | 1.47 | 1.46 | 1.45 |
FINAL STEP FOR FINDINGS: Making a Visual of Choice
Visual of choice
# visualize
Iceland_visual <- table(iceland_data_clean$Trust_in_Scientists, iceland_data_clean$Ideology) %>%
cprop() %>% ## this is the function to calculate column percentages
as.data.frame() %>%
filter(Var1 != "Total",
Var2 != "All") %>%
ggplot(aes(x=Var1, y=Freq, fill=Var2)) +
geom_col(position = "dodge") +
labs(title="Trust in Scientists in Iceland by Political Ideology",
y = "Conditional Percentage",
x = "No Trust at All to Complete Trust in Scientists",
fill = "Political Ideology") +
scale_fill_manual(values = c("Left" = "red",
"Right" = "blue",
"Centrist" = "lightblue"
))Displaying the visual
Iceland_visual
Saving the plot
ggsave(filename = "Iceland_visual", plot = Iceland_visual, device = "pdf", width = 10, height = 8)