Homework 3
Mingdi Zhou
2024-03-10
Setting up environment
# 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: 程辑包'infer'是用R版本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")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.
france_data <- ess %>%
# FiTter for France
filter(cntry == "FR") %>%
mutate(
# Recoding 'lrscale' into a new variable 'ID'
ID = case_when(
lrscale %in% 0:3 ~ "Left", # Values 0 to 3 in 'lrscale' are categorized as "Left"
lrscale %in% 7:10 ~ "Right", # Values 7 to 10 in 'lrscale' are categorized as "Right"
lrscale %in% c(4, 5, 6, 77, 88, 99) ~ NA_character_ # Values 77, 88, 99 in 'lrscale' are set as NA
),
# Recoding education with case_when for consistency
educ.ba = case_when(
essround < 5 & edulvla == 5 ~ "BA or more",
essround >= 5 & edulvlb > 600 ~ "BA or more",
TRUE ~ "No BA"
)
) %>%
# Ensure all NAs are handled uniformly
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) # Force educ.ba to be treated as character
)## Not necessary, but always good to rename when doing operations so you can backtrack to clean dataset
tab_dat <- france_data
# Create a crosstable and save it into the object mytab
mytab <- table(tab_dat$ID, tab_dat$educ.ba)
# Calculate row-wise sums of counts in table
rsums <- rowSums(mytab)
# Calculate column-wise sums of counts in table
csums <- colSums(mytab)
# Get the total count of the table (N)
N <- sum(mytab)
# Generate table of expected proportions
ptab <- tcrossprod(rsums/N, csums/N)
cat("Table of Expected Proportions:\n")## Table of Expected Proportions:print(ptab)## [,1] [,2]
## [1,] 0.1311243 0.3914592
## [2,] 0.1197912 0.3576253# Table of expected frequencies
ftab <- N * ptab
cat("Table of Expected Frequencies:\n")## Table of Expected Frequencies:print(ftab)## [,1] [,2]
## [1,] 1181.561 3527.439
## [2,] 1079.439 3222.561# significance level
alpha <- 0.05
# Critical Value for Independence:
c_val <- qchisq(1 - alpha, df = 1, lower.tail = FALSE)
cat("Critical Value:", round(c_val, 3), "\n")## Critical Value: 0.004# Pearson's chi-squared statistic:
test_stat <- sum((mytab - ftab)^2 / ftab)
cat("Pearson's X^2:", round(test_stat, 4), "\n")## Pearson's X^2: 100.8551# p-value:
p_val <- pchisq(test_stat, df = 1, lower.tail = FALSE)
cat("p-value:",round(p_val, 4), "\n")## p-value: 0# Using chisq.test for validation:
cat("Chi-squared test result:\n")## Chi-squared test result:## Chi-squared test result:
print(chisq.test(mytab, correct = F))##
## Pearson's Chi-squared test
##
## data: mytab
## X-squared = 100.86, df = 1, p-value < 2.2e-16Interpretation 1:
A Chi-Square statistic > 0.004 indicates a significant association between political orientation and education level at the 5% significance level, rejecting the null hypothesis of independence. If the statistic ≤ 0.004, the null hypothesis is not rejected, suggesting the observed relationship could be due to chance, lacking evidence of significance at the 5% level.
Interpretation 2:
The Chi-Square test with a statistic of 100.8551 reveals a significant discrepancy between observed frequencies and those expected if political orientation and education level were independent, indicating a strong statistical relationship. This high value suggests that distributions of education levels among political orientations deviate significantly from what would be expected by chance, implying an association between someone’s political orientation and their education level in the dataset.
Interpretation 3:
With a p-value of 0, there’s strong evidence of a significant association between political orientation and education level, suggesting the distribution of education varies with political leanings. This leads to rejecting the null hypothesis of independence, indicating a dependency between political orientation and the likelihood of having a BA degree.
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).
france_data <- france_data %>%
mutate(
sbbsntx = case_when(
sbbsntx == 1 ~ "Agree Strongly",
sbbsntx == 2 ~ "Agree",
sbbsntx == 3 ~ "Neither agree nor disagree",
sbbsntx == 4 ~ "Disagree",
sbbsntx == 5 ~ "Disagree strongly",
TRUE ~ NA_character_
),
domicil = case_when(
domicil == 1 ~ "Urban",
domicil == 2 ~ "Urban",
domicil == 3 ~ "Rural",
domicil == 4 ~ "Rural",
domicil == 5 ~ "Rural",
TRUE ~ NA_character_
)
)
france_data <- france_data %>%
filter( !is.na(sbbsntx),!is.na(domicil))test_stat <- france_data %>%
specify(explanatory = domicil, # change variable name for explanatory variable
response = sbbsntx) %>% # change variable name for outcome of interest
hypothesize(null = "independence") %>%
calculate(stat = "Chisq") # replace in between quotation marks appropriate test statistic
print(test_stat$stat)## X-squared
## 19.97187null_distribution <- france_data %>%
specify(explanatory = domicil,
response = sbbsntx) %>%
hypothesize(null = "independence") %>%
generate(reps = 1000, type = "permute") %>% # only line we are adding, use reps = 1000 as standard but change if too computationally demanding
calculate(stat = "Chisq")p_val <- null_distribution %>% # replace name here if you assigned something other than null_distribution above
get_pvalue(obs_stat = test_stat, direction = "two-sided") # would only replace test_stat if assigned another name in Step 1
p_val## # A tibble: 1 × 1
## p_value
## <dbl>
## 1 0.004Given the X-squared value of 19.97187 and P-value of 0.004, there’s a statistically significant association between domicile type (urban vs. rural) and sbbsntx responses in France
##Task 3
null_distribution %>%
visualize() +
shade_p_value(obs_stat = test_stat, direction = "two-sided")## 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 1.71
## 2 2 4.38
## 3 3 3.07
## 4 4 4.85
## 5 5 3.94
## 6 6 3.76
## 7 7 2.33
## 8 8 8.83
## 9 9 1.05
## 10 10 0.908
## # ℹ 990 more rowsconf_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 1.71
## 2 2 4.38
## 3 3 3.07
## 4 4 4.85
## 5 5 3.94
## 6 6 3.76
## 7 7 2.33
## 8 8 8.83
## 9 9 1.05
## 10 10 0.908
## # ℹ 990 more rowsInterpretation:
The graph is a right-skewed null distribution of chi-square test statistics from simulations under the null hypothesis of no association. The observed test statistic is about 19.97, indicated by the red line far into the right tail, suggesting a statistically significant result since it falls outside the common values. The associated p-value is 0.004, which is below the typical significance level of 0.05, leading to the rejection of the null hypothesis and the conclusion of a significant association between domicile type and survey responses.
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).
france_data <- france_data %>%
mutate(
ccrdprs = case_when(
ccrdprs %in% 0 ~ "Not at all",
ccrdprs %in% 1 ~ "1",
ccrdprs %in% 2 ~ "2",
ccrdprs %in% 3 ~ "3",
ccrdprs %in% 4 ~ "4",
ccrdprs %in% 5 ~ "5",
ccrdprs %in% 6 ~ "6",
ccrdprs %in% 7 ~ "7",
ccrdprs %in% 8 ~ "8",
ccrdprs %in% 9 ~ "9",
ccrdprs %in% 10 ~ "A great deal",
TRUE ~ NA_character_
),
lrscale = case_when(
lrscale %in% 0:3 ~ "Left", # Values 0 to 3 in 'lrscale' are categorized as "Left"
lrscale %in% 7:10 ~ "Right", # Values 7 to 10 in 'lrscale' are categorized as "Right"
lrscale %in% c(4, 5, 6, 77, 88, 99) ~ NA_character_ # Values 77, 88, 99 in 'lrscale' are set as NA
)
)
france_data<-france_data%>%
filter(!is.na(ccrdprs),!is.na(lrscale))test_stat <- france_data %>%
specify(explanatory = lrscale, # Ensure this is your explanatory variable
response = ccrdprs) %>% # Ensure this is your response variable
hypothesize(null = "independence") %>%
calculate(stat = "Chisq") # Chi-square statistic for testing independence
print(test_stat$stat)## X-squared
## 30.77499null_distribution <- france_data %>%
specify(explanatory = lrscale,
response = ccrdprs) %>%
hypothesize(null = "independence") %>%
generate(reps = 1000, type = "permute") %>% # only line we are adding, use reps = 1000 as standard but change if too computationally demanding
calculate(stat = "Chisq")p_val <- null_distribution %>% # replace name here if you assigned something other than null_distribution above
get_pvalue(obs_stat = test_stat, direction = "two-sided") # would only replace test_stat if assigned another name in Step 1## 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 0Task 5
null_distribution %>%
visualize() +
shade_p_value(obs_stat = test_stat, direction = "two-sided")## 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.7
## 2 2 15.9
## 3 3 11.7
## 4 4 5.14
## 5 5 15.6
## 6 6 4.96
## 7 7 14.3
## 8 8 3.76
## 9 9 5.93
## 10 10 8.52
## # ℹ 990 more rowsInterpretation:
The graph is a right-skewed null distribution of chi-square test statistics from simulations under the null hypothesis of no association. The observed test statistic is about 30.77, indicated by the red line far into the right tail, suggesting a statistically significant result since it falls outside the common values. The associated p-value is 0.002, which is below the typical significance level of 0.05, leading to the rejection of the null hypothesis and the conclusion of a significant association between lrscale type and survey responses.