Cappuccio_Julian_Homework_3
Julian Cappuccio
2024-03-11
# List of packages
packages <- c("tidyverse", "infer", "fst") # 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.5
## ✔ 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## Warning: package 'infer' was built under R version 4.3.3## [[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"# ess <- read_fst("All-ESS-Data.fst")# List of packages
packages <- c("tidyverse", "modelsummary", "flextable",
"fst", "viridis", "knitr", "kableExtra", "rmarkdown", "ggridges", "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 package: 'flextable'## The following object is masked from 'package:purrr':
##
## compose## Loading required package: viridisLite##
## 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] "fst" "infer" "lubridate" "forcats" "stringr" "dplyr"
## [7] "purrr" "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [13] "stats" "graphics" "grDevices" "utils" "datasets" "methods"
## [19] "base"
##
## [[2]]
## [1] "modelsummary" "fst" "infer" "lubridate" "forcats"
## [6] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [11] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [16] "grDevices" "utils" "datasets" "methods" "base"
##
## [[3]]
## [1] "flextable" "modelsummary" "fst" "infer" "lubridate"
## [6] "forcats" "stringr" "dplyr" "purrr" "readr"
## [11] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [16] "graphics" "grDevices" "utils" "datasets" "methods"
## [21] "base"
##
## [[4]]
## [1] "flextable" "modelsummary" "fst" "infer" "lubridate"
## [6] "forcats" "stringr" "dplyr" "purrr" "readr"
## [11] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [16] "graphics" "grDevices" "utils" "datasets" "methods"
## [21] "base"
##
## [[5]]
## [1] "viridis" "viridisLite" "flextable" "modelsummary" "fst"
## [6] "infer" "lubridate" "forcats" "stringr" "dplyr"
## [11] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [16] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [21] "datasets" "methods" "base"
##
## [[6]]
## [1] "knitr" "viridis" "viridisLite" "flextable" "modelsummary"
## [6] "fst" "infer" "lubridate" "forcats" "stringr"
## [11] "dplyr" "purrr" "readr" "tidyr" "tibble"
## [16] "ggplot2" "tidyverse" "stats" "graphics" "grDevices"
## [21] "utils" "datasets" "methods" "base"
##
## [[7]]
## [1] "kableExtra" "knitr" "viridis" "viridisLite" "flextable"
## [6] "modelsummary" "fst" "infer" "lubridate" "forcats"
## [11] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [16] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [21] "grDevices" "utils" "datasets" "methods" "base"
##
## [[8]]
## [1] "rmarkdown" "kableExtra" "knitr" "viridis" "viridisLite"
## [6] "flextable" "modelsummary" "fst" "infer" "lubridate"
## [11] "forcats" "stringr" "dplyr" "purrr" "readr"
## [16] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [21] "graphics" "grDevices" "utils" "datasets" "methods"
## [26] "base"
##
## [[9]]
## [1] "ggridges" "rmarkdown" "kableExtra" "knitr" "viridis"
## [6] "viridisLite" "flextable" "modelsummary" "fst" "infer"
## [11] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [16] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [21] "stats" "graphics" "grDevices" "utils" "datasets"
## [26] "methods" "base"
##
## [[10]]
## [1] "questionr" "ggridges" "rmarkdown" "kableExtra" "knitr"
## [6] "viridis" "viridisLite" "flextable" "modelsummary" "fst"
## [11] "infer" "lubridate" "forcats" "stringr" "dplyr"
## [16] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [21] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [26] "datasets" "methods" "base"ess <- read_fst("All-ESS-Data.fst")france_data <- read.fst("france_data.fst")Homework 3 (5%): due by next lecture on March 12th
Instructions: Start a new R markdown for the homework and call it “Yourlastname_Firstname_Homework_2”.
Copy everything below from Task 1 to Task 5. Keep the task prompt and questions, and provide your code and answer underneath.
Remember: you need all the steps for your code to work, including loading your data – otherwise it will not knit.
To generate a new code box, click on the +C sign above. Underneath your code, provide your answer to the task question.
When you are done, click on “Knit” above, then “Knit to Html”. Wait for everything to compile. If you get an error like “Execution halted”, it means there are issues with your code you must fix. When all issues are fixed, it will prompt a new window. Then click on “Publish” in the top right, and then Rpubs (the first option) and follow the instructions to create your Rpubs account and get your Rpubs link for your document (i.e., html link as I provide for the tutorial).
Note: Make sure to provide both your markdown file and R pubs link. If you do not submit both, you will be penalized 2 pts. out of the 5 pts. total.
Task 1
Based on the recoding of lrscale for “left” and “right” and omitting “moderates” (see Tutorial 4), and educ.ba from this tutorial, do the long form coding for Chi-Square. Using the steps outlined in the tutorial, first generate the tables of expected proportions and frequencies. Determine and interpret the critical value for independence. Finally, determine and interpret both the Pearson’s chi-squared statistic and the p-value. Discuss main takeaways.
Answer to Task 1
Recoding of lrscale and educ.ba
france_data <- france_data %>%
mutate(
lrscale = case_when(
lrscale %in% 0:3 ~ "Left",
lrscale %in% 7:10 ~ "Right", TRUE ~ NA_character_
),
)france_data <- ess %>%
filter(cntry == "FR") %>%
mutate(
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(55, 77, 88, 99), NA_integer_, edulvla),
edulvlb = ifelse(edulvlb %in% c(5555, 7777, 8888, 9999), NA_integer_, edulvlb),
educ.ba = as.character(educ.ba)
)Generate Table of proportions and frequencies
tab_dat <- france_data
mytab <- table(tab_dat$lrscale, tab_dat$educ.ba)
rsums <- rowSums(mytab)
csums <- colSums(mytab)
N <- sum(mytab)
ptab <- tcrossprod(rsums/N, csums/N)
cat("Table of expected proportions")## Table of expected proportionsprint(ptab)## [,1] [,2]
## [1,] 0.011086193 0.03855137
## [2,] 0.005853979 0.02035676
## [3,] 0.014887173 0.05176899
## [4,] 0.023415916 0.08142703
## [5,] 0.021562351 0.07498140
## [6,] 0.061871513 0.21515338
## [7,] 0.018699885 0.06502740
## [8,] 0.020518255 0.07135064
## [9,] 0.017104411 0.05947926
## [10,] 0.004798151 0.01668520
## [11,] 0.008047755 0.02798544
## [12,] 0.004469671 0.01554294
## [13,] 0.011027535 0.03834740ftab <- N * ptab
cat("Table of expected frequencies:\n")## Table of expected frequencies:print(ftab)## [,1] [,2]
## [1,] 211.0589 733.9411
## [2,] 111.4481 387.5519
## [3,] 283.4220 985.5780
## [4,] 445.7922 1550.2078
## [5,] 410.5040 1427.4960
## [6,] 1177.9099 4096.0901
## [7,] 356.0084 1237.9916
## [8,] 390.6265 1358.3735
## [9,] 325.6338 1132.3662
## [10,] 91.3472 317.6528
## [11,] 153.2132 532.7868
## [12,] 85.0936 295.9064
## [13,] 209.9422 730.0578alpha <- 0.05
c_val <- qchisq(alpha, df = 1, lower.tail = FALSE)
cat("Critical Value:", round(c_val, 3), "\n")## Critical Value: 3.841Based on this information, we can interpret that the critical value is 3.841.
Task 2
Do Steps 1 to 3 using the infer package (as we did in the tutorial), for the variables sbbsntx (recoded to the corresponding 5 categories) as the response and domicil (recoded to urban and rural, setting 2 to urban instead of peri-urban) as the explanatory variable for the country of France.
See variable info here: https://ess.sikt.no/en/datafile/ffc43f48-e15a-4a1c-8813-47eda377c355/92?tab=1&elements=[%22874a86e2-a0f6-40b4-aef7-51c8eed98a7d/1%22])
Provide interpretations of the output (consider the variable info).
Answer to Task 2
france_data <- ess %>%
filter(cntry == "FR") %>%
mutate(sbbsntx = case_when(sbbsntx == 1 ~ "Agree strongly", sbbsntx == 2 ~ "Agree", sbbsntx == 3 ~ "Neither agree nor disagree", sbbsntx == 4 ~ "Disagree", sbbsntx == 5 ~ "Disagree strongly", sbbsntx %in% c(7, 8, 9) ~ NA_character_, TRUE ~ as.character(sbbsntx)),
domicil = case_when(domicil == 1 ~ "Urban", domicil == 2 ~ "Urban", domicil == 3 ~ "Rural", domicil == 4 ~ "Rural", domicil == 5 ~ "Rural", domicil %in% c(7, 8, 9) ~ NA_character_, TRUE ~ as.character(domicil)))test_stat <- france_data %>%
specify(explanatory = domicil, response = sbbsntx) %>%
hypothesize(null = "independence") %>%
calculate(stat = "Chisq")## Warning: Removed 14999 rows containing missing values.print(test_stat$stat)## X-squared
## 19.97187null_distribution <- france_data %>%
specify(explanatory = domicil, response = sbbsntx) %>%
hypothesize(null = "independence") %>%
generate(reps = 1000, type = "permute") %>%
calculate(stat = "Chisq")## Warning: Removed 14999 rows containing missing values.p_val <- null_distribution %>% get_pvalue(obs_stat = test_stat, direction = "two-sided") ## 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.p_val## # A tibble: 1 × 1
## p_value
## <dbl>
## 1 0After completing the first 3 steps we are able to interpret the output. Working through these calculations has left us with a p value of 0. With this in mind, this means we would reject the null hypothesis as the P value for this question is below 0.05.
Task 3
Now do Step 4 (i.e., the null distribution visualization with confidence intervals) for the same variables and country as Task 2, and interpret the output.
Answer to Task 3
conf_int <- null_distribution %>%
get_confidence_interval(level = 0.95, type = "percentile")
null_distribution %>%
visualize() + shade_p_value(obs_stat = test_stat, direction = "two-sided") + shade_confidence_interval(endpoints = conf_int)## Warning: Chi-Square usually corresponds to right-tailed tests. Proceed with
## caution.
null_distribution## Response: sbbsntx (factor)
## Explanatory: domicil (factor)
## Null Hypothesis: independence
## # A tibble: 1,000 × 2
## replicate stat
## <int> <dbl>
## 1 1 7.25
## 2 2 0.641
## 3 3 3.71
## 4 4 4.61
## 5 5 1.76
## 6 6 2.57
## 7 7 2.01
## 8 8 5.69
## 9 9 11.5
## 10 10 1.09
## # ℹ 990 more rowsWith this visualization, we are now able to interpret this output. This visualization shows and confirms that we should reject the null hypothesis.
Task 4
Conduct Steps 1 to 3 using the infer package (as we did in the tutorial), for the variables ccrdprs (leaving it on the 0-10 numeric scale) as the response (or outcome) variable and lrscale (recoded as left and right, omitting “moderates” as we did in Task 1) as the explanatory variable.
Variable info here: https://ess.sikt.no/en/datafile/ffc43f48-e15a-4a1c-8813-47eda377c355/92?tab=1&elements=[%2283b08cf2-508e-49a0-9fc8-c3cc281290ec/4%22]
Provide interpretations of the output (consider the variable info).
Answer to Task 4
france_data <- ess %>%
filter(cntry == "FR") %>%
mutate(ccrdprs = case_when(ccrdprs %in% c(66, 77, 88, 99) ~ NA_character_, TRUE ~ as.character(ccrdprs)),
lrscale = case_when(lrscale %in% 0:3 ~ "Left", lrscale %in% 7:10 ~ "Right", TRUE ~ NA_character_
)
)test_stat <- france_data %>%
specify(explanatory = lrscale, response = ccrdprs) %>%
hypothesize(null = "independence") %>%
calculate(stat = "Chisq")## Warning: Removed 17167 rows containing missing values.print(test_stat$stat)## X-squared
## 35.04425null_distribution <- france_data %>%
specify(explanatory = lrscale, response = ccrdprs) %>%
hypothesize(null = "independence") %>%
generate(reps = 1000, type = "permute") %>%
calculate(stat = "Chisq")## Warning: Removed 17167 rows containing missing values.p_val <- null_distribution %>% get_pvalue(obs_stat = test_stat, direction = "two-sided") ## 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.p_val## # A tibble: 1 × 1
## p_value
## <dbl>
## 1 0Based on this information, the p value that we get is a p value of 0. This means that the null hypothesis can be rejected.
Task 5
Now do Step 4 (i.e., the null distribution visualization with confidence intervals) for the same variables and country as Task 4, and interpret the output.
conf_int <- null_distribution %>%
get_confidence_interval(level = 0.95, type = "percentile")
null_distribution %>%
visualize() + shade_p_value(obs_stat = test_stat, direction = "two-sided") + shade_confidence_interval(endpoints = conf_int)## Warning: Chi-Square usually corresponds to right-tailed tests. Proceed with
## caution.
null_distribution## Response: ccrdprs (factor)
## Explanatory: lrscale (factor)
## Null Hypothesis: independence
## # A tibble: 1,000 × 2
## replicate stat
## <int> <dbl>
## 1 1 14.2
## 2 2 5.85
## 3 3 6.88
## 4 4 11.3
## 5 5 10.8
## 6 6 9.80
## 7 7 7.46
## 8 8 5.63
## 9 9 12.1
## 10 10 6.76
## # ℹ 990 more rowsBased on this visualization, we will reject the null.