setwd("C:/Users/Owner/Downloads/All Downloads/Stepping 2 Workspace")
# List of packages
packages <- c("tidyverse", "modelsummary", "forcats", "RColorBrewer",
"fst", "viridis", "knitr", "rmarkdown", "ggridges", "viridis", "questionr", "flextable", "infer", "broom", "htmltools") # 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: viridisLite
##
##
## Attaching package: 'flextable'
##
##
## The following object is masked from 'package:purrr':
##
## compose
## Warning: package 'htmltools' was built under R version 4.3.2
## [[1]]
## [1] "lubridate" "forcats" "stringr" "dplyr" "purrr" "readr"
## [7] "tidyr" "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [13] "grDevices" "utils" "datasets" "methods" "base"
##
## [[2]]
## [1] "modelsummary" "lubridate" "forcats" "stringr" "dplyr"
## [6] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [11] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [16] "datasets" "methods" "base"
##
## [[3]]
## [1] "modelsummary" "lubridate" "forcats" "stringr" "dplyr"
## [6] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [11] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [16] "datasets" "methods" "base"
##
## [[4]]
## [1] "RColorBrewer" "modelsummary" "lubridate" "forcats" "stringr"
## [6] "dplyr" "purrr" "readr" "tidyr" "tibble"
## [11] "ggplot2" "tidyverse" "stats" "graphics" "grDevices"
## [16] "utils" "datasets" "methods" "base"
##
## [[5]]
## [1] "fst" "RColorBrewer" "modelsummary" "lubridate" "forcats"
## [6] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [11] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [16] "grDevices" "utils" "datasets" "methods" "base"
##
## [[6]]
## [1] "viridis" "viridisLite" "fst" "RColorBrewer" "modelsummary"
## [6] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [11] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [16] "stats" "graphics" "grDevices" "utils" "datasets"
## [21] "methods" "base"
##
## [[7]]
## [1] "knitr" "viridis" "viridisLite" "fst" "RColorBrewer"
## [6] "modelsummary" "lubridate" "forcats" "stringr" "dplyr"
## [11] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [16] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [21] "datasets" "methods" "base"
##
## [[8]]
## [1] "rmarkdown" "knitr" "viridis" "viridisLite" "fst"
## [6] "RColorBrewer" "modelsummary" "lubridate" "forcats" "stringr"
## [11] "dplyr" "purrr" "readr" "tidyr" "tibble"
## [16] "ggplot2" "tidyverse" "stats" "graphics" "grDevices"
## [21] "utils" "datasets" "methods" "base"
##
## [[9]]
## [1] "ggridges" "rmarkdown" "knitr" "viridis" "viridisLite"
## [6] "fst" "RColorBrewer" "modelsummary" "lubridate" "forcats"
## [11] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [16] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [21] "grDevices" "utils" "datasets" "methods" "base"
##
## [[10]]
## [1] "ggridges" "rmarkdown" "knitr" "viridis" "viridisLite"
## [6] "fst" "RColorBrewer" "modelsummary" "lubridate" "forcats"
## [11] "stringr" "dplyr" "purrr" "readr" "tidyr"
## [16] "tibble" "ggplot2" "tidyverse" "stats" "graphics"
## [21] "grDevices" "utils" "datasets" "methods" "base"
##
## [[11]]
## [1] "questionr" "ggridges" "rmarkdown" "knitr" "viridis"
## [6] "viridisLite" "fst" "RColorBrewer" "modelsummary" "lubridate"
## [11] "forcats" "stringr" "dplyr" "purrr" "readr"
## [16] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [21] "graphics" "grDevices" "utils" "datasets" "methods"
## [26] "base"
##
## [[12]]
## [1] "flextable" "questionr" "ggridges" "rmarkdown" "knitr"
## [6] "viridis" "viridisLite" "fst" "RColorBrewer" "modelsummary"
## [11] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [16] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [21] "stats" "graphics" "grDevices" "utils" "datasets"
## [26] "methods" "base"
##
## [[13]]
## [1] "infer" "flextable" "questionr" "ggridges" "rmarkdown"
## [6] "knitr" "viridis" "viridisLite" "fst" "RColorBrewer"
## [11] "modelsummary" "lubridate" "forcats" "stringr" "dplyr"
## [16] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [21] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [26] "datasets" "methods" "base"
##
## [[14]]
## [1] "broom" "infer" "flextable" "questionr" "ggridges"
## [6] "rmarkdown" "knitr" "viridis" "viridisLite" "fst"
## [11] "RColorBrewer" "modelsummary" "lubridate" "forcats" "stringr"
## [16] "dplyr" "purrr" "readr" "tidyr" "tibble"
## [21] "ggplot2" "tidyverse" "stats" "graphics" "grDevices"
## [26] "utils" "datasets" "methods" "base"
##
## [[15]]
## [1] "htmltools" "broom" "infer" "flextable" "questionr"
## [6] "ggridges" "rmarkdown" "knitr" "viridis" "viridisLite"
## [11] "fst" "RColorBrewer" "modelsummary" "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")
# Filter for France
# Recoding outcome variables, response= trstplc, explanatory= gender= gndr, edulvla, paid work status=pdwrk
finland_data <- ess %>% filter(cntry == "FI") %>% mutate(trstplc = ifelse(trstplc %in% c(77, 88, 99), NA, trstplc),
# Recoding gender
gndr = case_when(
gndr == 1 ~ "Male",
gndr == 2 ~ "Female",
gndr == 9 ~ NA_character_,
TRUE ~ as.character(gndr)
),
# Recoding paid work status
pdwrk = case_when(
pdwrk == 1 ~ "Yes",
pdwrk == 0 ~ "No",
pdwrk == 9 ~ NA_character_,
TRUE ~ as.character(gndr)
)
)
#Plotting outcome and explanatory.
ggplot(finland_data, mapping = aes(x = gndr, fill = trstplc)) + geom_bar()
## Warning: The following aesthetics were dropped during statistical transformation: fill
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
## the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
## variable into a factor?
#Hypothesize null of outcome and explanatory
finland_data_shuffled <- finland_data %>%
group_by(gndr, trstplc) %>%
tally()
ggplot(finland_data_shuffled, mapping = aes(x = gndr, fill = trstplc)) + geom_bar()
## Warning: The following aesthetics were dropped during statistical transformation: fill
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
## the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
## variable into a factor?
# using unique to identify
unique(finland_data$gndr)
## [1] "Female" "Male"
# does not return NA, but the following does
unique(finland_data$trstplc)
## [1] 9 7 8 6 10 4 NA 5 3 1 0 2
# here's how to make sure it's really cleaned
finland_data <- finland_data %>% filter(!is.na(trstplc))
# double-checking + quick distribution
unique(finland_data$trstplc)
## [1] 9 7 8 6 10 4 5 3 1 0 2
table(finland_data$trstplc)
##
## 0 1 2 3 4 5 6 7 8 9 10
## 100 75 120 187 330 738 942 2325 5474 6554 2643
# Step 1: Calculate the test statistic of your sample
test_stat <- finland_data %>%
specify(explanatory = gndr,
response = trstplt) %>%
hypothesize(null = "independence") %>%
calculate(stat = "t")
## Warning: The statistic is based on a difference or ratio; by default, for
## difference-based statistics, the explanatory variable is subtracted in the
## order "Female" - "Male", or divided in the order "Female" / "Male" for
## ratio-based statistics. To specify this order yourself, supply `order =
## c("Female", "Male")` to the calculate() function.
print(test_stat$stat)
## t
## 5.185316
# Step 2: Simulate the null distribution
#Null hypothesis= Men and Women have the same rate of trust in the police
#Alternative hypothesis= Women does not have equal rate of trust in the police than man. H0 not equal 0, H0 > 0, H0 < 0.
null_distribution_gender<- finland_data %>%
specify(explanatory = gndr,
response = trstplc) %>%
hypothesize(null = "independence") %>%
generate(reps = 1000, type = "permute") %>%
calculate( stat = "t", order = c("Male", "Female"))
null_distribution_gender
## Response: trstplc (numeric)
## Explanatory: gndr (factor)
## Null Hypothesis: independence
## # A tibble: 1,000 × 2
## replicate stat
## <int> <dbl>
## 1 1 0.151
## 2 2 -0.536
## 3 3 0.813
## 4 4 0.727
## 5 5 1.17
## 6 6 0.0912
## 7 7 0.0740
## 8 8 0.263
## 9 9 0.581
## 10 10 0.444
## # ℹ 990 more rows
# Step 3: Calculate the p-value of your sample
p_val_gender <- null_distribution_gender %>%
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()` for more information.
p_val_gender
## # A tibble: 1 × 1
## p_value
## <dbl>
## 1 0
# Step 4: Visualize
#Vizualizing confidence intervals
conf_int <- null_distribution_gender %>%
get_confidence_interval(level = 0.95, type = "percentile")
null_distribution_gender %>%
visualize() +
shade_p_value(obs_stat = test_stat, direction = "two-sided") +
shade_confidence_interval(endpoints = conf_int)
null_distribution_gender
## Response: trstplc (numeric)
## Explanatory: gndr (factor)
## Null Hypothesis: independence
## # A tibble: 1,000 × 2
## replicate stat
## <int> <dbl>
## 1 1 0.151
## 2 2 -0.536
## 3 3 0.813
## 4 4 0.727
## 5 5 1.17
## 6 6 0.0912
## 7 7 0.0740
## 8 8 0.263
## 9 9 0.581
## 10 10 0.444
## # ℹ 990 more rows
#Modelling with LM because trust in police is not binary, and this is the line of code in R to produce a linear regression model.
model1 <- lm(trstplc ~ gndr, data = finland_data)
# The coefficients
coefficients1 <- coef(model1)
print(coefficients1)
## (Intercept) gndrMale
## 8.08618619 -0.05923314
#My own interpretation
#(intercept of education level)
#if a person is a male= 8.08618619
#If a person is a female= 8.08618619 -0.05923314 6
#The intercept represents the predicted value of the outcome variable (in this case, Trust in Police 0-10) when the explanatory variable is zero or, in the case of categorical variables, at their reference level.
#In simpler terms, when looking at the category that's "Male" in Gender (thus, holding a "Female"), the expected value of trstplc or trust in police is approximately 8.08618619 -0.05923314 = .
# Table presentation
tidy(model1) |>
knitr::kable (digits = 3)
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | 8.086 | 0.017 | 484.919 | 0.000 |
gndrMale | -0.059 | 0.024 | -2.480 | 0.013 |
## Three models
#The model1 is a predictor, model2 represents a second predictor, and model3 is an incorporation of an interaction.
# paid work status co variate or predictor is suggested by the interaction mode because the first model is one predictor, the second is two predictors, and third incorporates an interaction. The idea is female does not have equal rates of trust in police than male, which is further explored by paid work status (Yes to having Paid Work Status vs No to having Paid Work Status).Looking at the coefficients, standard errors, and adjust R-squad is importat for interpretations.
# List of packages
packages <- c("tidyverse", "infer", "fst", "modelsummary", "effects", "survey", "MASS") # 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)
## 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: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
## [[1]]
## [1] "fstcore" "htmltools" "broom" "infer" "flextable"
## [6] "questionr" "ggridges" "rmarkdown" "knitr" "viridis"
## [11] "viridisLite" "fst" "RColorBrewer" "modelsummary" "lubridate"
## [16] "forcats" "stringr" "dplyr" "purrr" "readr"
## [21] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [26] "graphics" "grDevices" "utils" "datasets" "methods"
## [31] "base"
##
## [[2]]
## [1] "fstcore" "htmltools" "broom" "infer" "flextable"
## [6] "questionr" "ggridges" "rmarkdown" "knitr" "viridis"
## [11] "viridisLite" "fst" "RColorBrewer" "modelsummary" "lubridate"
## [16] "forcats" "stringr" "dplyr" "purrr" "readr"
## [21] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [26] "graphics" "grDevices" "utils" "datasets" "methods"
## [31] "base"
##
## [[3]]
## [1] "fstcore" "htmltools" "broom" "infer" "flextable"
## [6] "questionr" "ggridges" "rmarkdown" "knitr" "viridis"
## [11] "viridisLite" "fst" "RColorBrewer" "modelsummary" "lubridate"
## [16] "forcats" "stringr" "dplyr" "purrr" "readr"
## [21] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [26] "graphics" "grDevices" "utils" "datasets" "methods"
## [31] "base"
##
## [[4]]
## [1] "fstcore" "htmltools" "broom" "infer" "flextable"
## [6] "questionr" "ggridges" "rmarkdown" "knitr" "viridis"
## [11] "viridisLite" "fst" "RColorBrewer" "modelsummary" "lubridate"
## [16] "forcats" "stringr" "dplyr" "purrr" "readr"
## [21] "tidyr" "tibble" "ggplot2" "tidyverse" "stats"
## [26] "graphics" "grDevices" "utils" "datasets" "methods"
## [31] "base"
##
## [[5]]
## [1] "effects" "carData" "fstcore" "htmltools" "broom"
## [6] "infer" "flextable" "questionr" "ggridges" "rmarkdown"
## [11] "knitr" "viridis" "viridisLite" "fst" "RColorBrewer"
## [16] "modelsummary" "lubridate" "forcats" "stringr" "dplyr"
## [21] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [26] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [31] "datasets" "methods" "base"
##
## [[6]]
## [1] "survey" "survival" "Matrix" "grid" "effects"
## [6] "carData" "fstcore" "htmltools" "broom" "infer"
## [11] "flextable" "questionr" "ggridges" "rmarkdown" "knitr"
## [16] "viridis" "viridisLite" "fst" "RColorBrewer" "modelsummary"
## [21] "lubridate" "forcats" "stringr" "dplyr" "purrr"
## [26] "readr" "tidyr" "tibble" "ggplot2" "tidyverse"
## [31] "stats" "graphics" "grDevices" "utils" "datasets"
## [36] "methods" "base"
##
## [[7]]
## [1] "MASS" "survey" "survival" "Matrix" "grid"
## [6] "effects" "carData" "fstcore" "htmltools" "broom"
## [11] "infer" "flextable" "questionr" "ggridges" "rmarkdown"
## [16] "knitr" "viridis" "viridisLite" "fst" "RColorBrewer"
## [21] "modelsummary" "lubridate" "forcats" "stringr" "dplyr"
## [26] "purrr" "readr" "tidyr" "tibble" "ggplot2"
## [31] "tidyverse" "stats" "graphics" "grDevices" "utils"
## [36] "datasets" "methods" "base"
model1 <- lm(trstplc ~ gndr, data = finland_data)
model2 <- lm(trstplc ~ gndr + pdwrk, data = finland_data)
model3 <- lm(trstplc ~ gndr + pdwrk + gndr*pdwrk, data = finland_data)
modelsummary(
list(model1, model2, model3))
(1) | (2) | (3) | |
---|---|---|---|
(Intercept) | 8.086 | 8.014 | 8.042 |
(0.017) | (0.020) | (0.023) | |
gndrMale | −0.059 | −0.070 | −0.132 |
(0.024) | (0.024) | (0.035) | |
pdwrkYes | 0.148 | 0.091 | |
(0.024) | (0.033) | ||
gndrMale × pdwrkYes | 0.118 | ||
(0.048) | |||
Num.Obs. | 19488 | 19488 | 19488 |
R2 | 0.000 | 0.002 | 0.003 |
R2 Adj. | 0.000 | 0.002 | 0.002 |
AIC | 75219.2 | 75183.1 | 75179.0 |
BIC | 75242.9 | 75214.6 | 75218.4 |
Log.Lik. | −37606.625 | −37587.537 | −37584.517 |
RMSE | 1.67 | 1.66 | 1.66 |
# Now, let's extract just the coefficients
coefficients1 <- coef(model1)
print(coefficients1)
## (Intercept) gndrMale
## 8.08618619 -0.05923314
coefficients2 <- coef(model2)
print(coefficients2)
## (Intercept) gndrMale pdwrkYes
## 8.01435519 -0.07037206 0.14801807
coefficients3 <- coef(model3)
print(coefficients3)
## (Intercept) gndrMale pdwrkYes gndrMale:pdwrkYes
## 8.04200700 -0.13208846 0.09103755 0.11775164
# Striving to improve table presentation
tidy(model1) |>
knitr::kable (digits = 3)
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | 8.086 | 0.017 | 484.919 | 0.000 |
gndrMale | -0.059 | 0.024 | -2.480 | 0.013 |
tidy(model2) |>
knitr::kable(digits = 3)
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | 8.014 | 0.020 | 394.560 | 0.000 |
gndrMale | -0.070 | 0.024 | -2.941 | 0.003 |
pdwrkYes | 0.148 | 0.024 | 6.181 | 0.000 |
tidy(model3) |>
knitr::kable(digits = 3)
term | estimate | std.error | statistic | p.value |
---|---|---|---|---|
(Intercept) | 8.042 | 0.023 | 346.373 | 0.000 |
gndrMale | -0.132 | 0.035 | -3.808 | 0.000 |
pdwrkYes | 0.091 | 0.033 | 2.731 | 0.006 |
gndrMale:pdwrkYes | 0.118 | 0.048 | 2.458 | 0.014 |
# Using modelsummary for prettier table, focusing here on just R squared and adj.
modelsummary::modelsummary(
list(model1, model2, model3),
gof_map = c("nobs", "r.squared", "adj.r.squared"))
(1) | (2) | (3) | |
---|---|---|---|
(Intercept) | 8.086 | 8.014 | 8.042 |
(0.017) | (0.020) | (0.023) | |
gndrMale | −0.059 | −0.070 | −0.132 |
(0.024) | (0.024) | (0.035) | |
pdwrkYes | 0.148 | 0.091 | |
(0.024) | (0.033) | ||
gndrMale × pdwrkYes | 0.118 | ||
(0.048) | |||
Num.Obs. | 19488 | 19488 | 19488 |
R2 | 0.000 | 0.002 | 0.003 |
R2 Adj. | 0.000 | 0.002 | 0.002 |
#My own interpretation of Model2, the second predictor.
#(intercept of education level)
#if a person is a male= 8.014
#If a person is a female= 8.014 4 −0.070
#If a person has paid work= 8.014
#If a person has no paid work= 8.014 + 0.148
#standard error of model 1 (0.017). model 2 (0.020), and model 3(0.023) is associated with the intercept because this determines how uncertainty we can be around estimated intercept value.
#standard error of model 1 (0.024), model 2(0.024), and model 3(0.035) is associated with "gndrMale" while this determines certainty around the estimate difference.
#The numbers of observations being 19488 represents the number of observations used in each model, in other words, the sample size.
#In model 1, the R-squared (R^2) being the value of 0 represents the proportion of the variance in the outcome variable of the trust in police that's explained by the explanatory variable in the model.
#A Value r2= 0.002 or 0.003 means that only 0.2% or 0.3% of the variance in trust in police is explained by the paid work status. Here, the r2 value is considered very tiny ammount suggesting the paid work status alone would not explain much of the variation in trust in police.
#R2 Adj. 0.002 is the adjust r-squad meaning a modded version of R2 that utlizes the number of predictors from the model such that it is useful during the comparison of models from different numbers of predictors and they are similar to R2, which make sense as the model represent single predictor of paid work status.
#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) | |
---|---|---|---|
gndrMale | −0.1 (0.0)* | −0.1 (0.0)** | −0.1 (0.0)*** |
pdwrkYes | 0.1 (0.0)*** | 0.1 (0.0)** | |
gndrMale × pdwrkYes | 0.1 (0.0)* | ||
Num.Obs. | 19488 | 19488 | 19488 |
R2 | 0.000 | 0.002 | 0.003 |
R2 Adj. | 0.000 | 0.002 | 0.002 |
AIC | 75219.2 | 75183.1 | 75179.0 |
BIC | 75242.9 | 75214.6 | 75218.4 |
Log.Lik. | −37606.625 | −37587.537 | −37584.517 |
RMSE | 1.67 | 1.66 | 1.66 |
#The model summary belows represents what I have included in the model showing 1. gender of male, 2. having Paid Work Status, and 3. having both paid work status and gender of male.
### Seperate here
# We use the "effect" function to call the interaction term from model 1.
interaction_plot_1 <- effect("gndr", model1, na.rm=TRUE)
#We plot here the interaction term 1
plot(interaction_plot_1,
main="Interaction effect Female and Male Gender x Paid Work Status", ylab="Trust in Police",
xlab="Male and Female Gender")
# We use the "effect" function to call the interaction term from model 2.
interaction_plot_2 <- effect("gndr*pdwrk", model2, na.rm=TRUE)
## NOTE: gndr:pdwrk does not appear in the model
#We plot here the interaction term 2
plot(interaction_plot_2,
main="Interaction effect Female and Male Gender x Paid Work Status", ylab="Trust in Police",
xlab="Male and Female Gender")
# We use the "effect" function to call the interaction term from model 3.
interaction_plot_3 <- effect("gndr*pdwrk", model3, na.rm=TRUE)
#We plot here the interaction term 3
plot(interaction_plot_3,
main="Interaction effect Female and Male Gender x Paid Work Status", ylab="Trust in Police",
xlab="Male and Female Gender")
# Based on the interaction effect of gender(female and male) with trust in police, Model 3 zooms in the paid work status make it easier to interpret through a visualization compared to mode 1 or 2 which all shows the interaction effect of having paid work status for male and female does indeed have a positive effect on the trust in police.There is also difference of trust in police people depending on having paid work status, those that say yes to having paid work status has a tendency to trust in police more than those who say no to having paid work status.However, the relationship between gender and trust in police is already following a similar trend,the visualization indicates a tendency for Female to trust in police more than male, so regardless of having paid work status than it should not have an impact on the this gender and trust relationship. AS such, the models of the gender of female and male not having equal trust in police fits in with the alternative, which suggest evidence that might disprove the null hypothesis that Female and Male both have the rate of trust in police.
## Finally, Getting equation output with equatiomatic
remotes::install_github("datalorax/equatiomatic", force = TRUE)
## Downloading GitHub repo datalorax/equatiomatic@HEAD
## xfun (0.40 -> 0.41 ) [CRAN]
## rlang (1.1.1 -> 1.1.2 ) [CRAN]
## vctrs (0.6.3 -> 0.6.4 ) [CRAN]
## utf8 (1.2.3 -> 1.2.4 ) [CRAN]
## fansi (1.0.4 -> 1.0.5 ) [CRAN]
## stringi (1.7.12 -> 1.8.2 ) [CRAN]
## stringr (1.5.0 -> 1.5.1 ) [CRAN]
## dplyr (1.1.3 -> 1.1.4 ) [CRAN]
## httpuv (1.6.11 -> 1.6.12) [CRAN]
## knitr (1.44 -> 1.45 ) [CRAN]
## Installing 10 packages: xfun, rlang, vctrs, utf8, fansi, stringi, stringr, dplyr, httpuv, knitr
## Warning: packages 'stringr', 'dplyr', 'knitr' are in use and will not be
## installed
## Installing packages into 'C:/Users/Owner/AppData/Local/R/win-library/4.3'
## (as 'lib' is unspecified)
## package 'xfun' successfully unpacked and MD5 sums checked
## Warning: cannot remove prior installation of package 'xfun'
## Warning in file.copy(savedcopy, lib, recursive = TRUE): problem copying
## C:\Users\Owner\AppData\Local\R\win-library\4.3\00LOCK\xfun\libs\x64\xfun.dll to
## C:\Users\Owner\AppData\Local\R\win-library\4.3\xfun\libs\x64\xfun.dll:
## Permission denied
## Warning: restored 'xfun'
## package 'rlang' successfully unpacked and MD5 sums checked
## Warning: cannot remove prior installation of package 'rlang'
## Warning in file.copy(savedcopy, lib, recursive = TRUE): problem copying
## C:\Users\Owner\AppData\Local\R\win-library\4.3\00LOCK\rlang\libs\x64\rlang.dll
## to C:\Users\Owner\AppData\Local\R\win-library\4.3\rlang\libs\x64\rlang.dll:
## Permission denied
## Warning: restored 'rlang'
## package 'vctrs' successfully unpacked and MD5 sums checked
## Warning: cannot remove prior installation of package 'vctrs'
## Warning in file.copy(savedcopy, lib, recursive = TRUE): problem copying
## C:\Users\Owner\AppData\Local\R\win-library\4.3\00LOCK\vctrs\libs\x64\vctrs.dll
## to C:\Users\Owner\AppData\Local\R\win-library\4.3\vctrs\libs\x64\vctrs.dll:
## Permission denied
## Warning: restored 'vctrs'
## package 'utf8' successfully unpacked and MD5 sums checked
## Warning: cannot remove prior installation of package 'utf8'
## Warning in file.copy(savedcopy, lib, recursive = TRUE): problem copying
## C:\Users\Owner\AppData\Local\R\win-library\4.3\00LOCK\utf8\libs\x64\utf8.dll to
## C:\Users\Owner\AppData\Local\R\win-library\4.3\utf8\libs\x64\utf8.dll:
## Permission denied
## Warning: restored 'utf8'
## package 'fansi' successfully unpacked and MD5 sums checked
## Warning: cannot remove prior installation of package 'fansi'
## Warning in file.copy(savedcopy, lib, recursive = TRUE): problem copying
## C:\Users\Owner\AppData\Local\R\win-library\4.3\00LOCK\fansi\libs\x64\fansi.dll
## to C:\Users\Owner\AppData\Local\R\win-library\4.3\fansi\libs\x64\fansi.dll:
## Permission denied
## Warning: restored 'fansi'
## package 'stringi' successfully unpacked and MD5 sums checked
## Warning: cannot remove prior installation of package 'stringi'
## Warning in file.copy(savedcopy, lib, recursive = TRUE): problem copying
## C:\Users\Owner\AppData\Local\R\win-library\4.3\00LOCK\stringi\libs\icudt69l.dat
## to C:\Users\Owner\AppData\Local\R\win-library\4.3\stringi\libs\icudt69l.dat:
## Invalid argument
## Warning in file.copy(savedcopy, lib, recursive = TRUE): problem copying
## C:\Users\Owner\AppData\Local\R\win-library\4.3\00LOCK\stringi\libs\x64\stringi.dll
## to C:\Users\Owner\AppData\Local\R\win-library\4.3\stringi\libs\x64\stringi.dll:
## Permission denied
## Warning: restored 'stringi'
## package 'httpuv' successfully unpacked and MD5 sums checked
## Warning: cannot remove prior installation of package 'httpuv'
## Warning in file.copy(savedcopy, lib, recursive = TRUE): problem copying
## C:\Users\Owner\AppData\Local\R\win-library\4.3\00LOCK\httpuv\libs\x64\httpuv.dll
## to C:\Users\Owner\AppData\Local\R\win-library\4.3\httpuv\libs\x64\httpuv.dll:
## Permission denied
## Warning: restored 'httpuv'
##
## The downloaded binary packages are in
## C:\Users\Owner\AppData\Local\Temp\RtmpYn5q88\downloaded_packages
## ── R CMD build ─────────────────────────────────────────────────────────────────
##
checking for file 'C:\Users\Owner\AppData\Local\Temp\RtmpYn5q88\remotes2e3c3d3c4232\datalorax-equatiomatic-29ff168/DESCRIPTION' ...
checking for file 'C:\Users\Owner\AppData\Local\Temp\RtmpYn5q88\remotes2e3c3d3c4232\datalorax-equatiomatic-29ff168/DESCRIPTION' ...
✔ checking for file 'C:\Users\Owner\AppData\Local\Temp\RtmpYn5q88\remotes2e3c3d3c4232\datalorax-equatiomatic-29ff168/DESCRIPTION' (501ms)
##
─ preparing 'equatiomatic':
## checking DESCRIPTION meta-information ...
checking DESCRIPTION meta-information ...
✔ checking DESCRIPTION meta-information
##
─ checking for LF line-endings in source and make files and shell scripts
##
─ checking for empty or unneeded directories
##
─ building 'equatiomatic_0.3.1.tar.gz'
##
##
## Installing package into 'C:/Users/Owner/AppData/Local/R/win-library/4.3'
## (as 'lib' is unspecified)
equatiomatic::extract_eq(model1, use_coefs = TRUE)
\[ \operatorname{\widehat{trstplc}} = 8.09 - 0.06(\operatorname{gndr}_{\operatorname{Male}}) \]
equatiomatic::extract_eq(model2, use_coefs = TRUE)
\[ \operatorname{\widehat{trstplc}} = 8.01 - 0.07(\operatorname{gndr}_{\operatorname{Male}}) + 0.15(\operatorname{pdwrk}_{\operatorname{Yes}}) \]
equatiomatic::extract_eq(model3, use_coefs = TRUE)
\[ \operatorname{\widehat{trstplc}} = 8.04 - 0.13(\operatorname{gndr}_{\operatorname{Male}}) + 0.09(\operatorname{pdwrk}_{\operatorname{Yes}}) + 0.12(\operatorname{gndr}_{\operatorname{Male}} \times \operatorname{pdwrk}_{\operatorname{Yes}}) \]
Important note: look at the html link, it produces a nice written out version of the equation (the intended version) that is more interpretable than the output from just the markdown.
Please ask questions about your projects!
#copy and paste equation to quicklatex.com for the regression equation
Note that the echo = FALSE
parameter was added to the
code chunk to prevent printing of the R code that generated the
plot.