Nawaz_Fatima_Homework

Loading Packages

# List of packages
packages <- c("tidyverse", "infer", "fst", "modelsummary", "effects", "survey", "MASS", "aod", "interactions", "kableExtra", "flextable", "scales", "ggeffects") # 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.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: 'MASS'
## 
## 
## The following object is masked from 'package:dplyr':
## 
##     select
## 
## 
## 
## Attaching package: 'aod'
## 
## 
## The following object is masked from 'package:survival':
## 
##     rats
## 
## 
## 
## Attaching package: 'kableExtra'
## 
## 
## The following object is masked from 'package:dplyr':
## 
##     group_rows
## 
## 
## 
## Attaching package: 'flextable'
## 
## 
## The following objects are masked from 'package:kableExtra':
## 
##     as_image, footnote
## 
## 
## 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
## 
## 
## 
## Attaching package: 'ggeffects'
## 
## 
## The following object is masked from 'package:interactions':
## 
##     johnson_neyman

## [[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] "MASS"         "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] "aod"          "MASS"         "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] "interactions" "aod"          "MASS"         "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] "kableExtra"   "interactions" "aod"          "MASS"         "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"        
## 
## [[11]]
##  [1] "flextable"    "kableExtra"   "interactions" "aod"          "MASS"        
##  [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"        
## 
## [[12]]
##  [1] "scales"       "flextable"    "kableExtra"   "interactions" "aod"         
##  [6] "MASS"         "survey"       "survival"     "Matrix"       "grid"        
## [11] "effects"      "carData"      "modelsummary" "fst"          "infer"       
## [16] "lubridate"    "forcats"      "stringr"      "dplyr"        "purrr"       
## [21] "readr"        "tidyr"        "tibble"       "ggplot2"      "tidyverse"   
## [26] "stats"        "graphics"     "grDevices"    "utils"        "datasets"    
## [31] "methods"      "base"        
## 
## [[13]]
##  [1] "ggeffects"    "scales"       "flextable"    "kableExtra"   "interactions"
##  [6] "aod"          "MASS"         "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"

Task 1

Use data for Germany and the variable coding for model3 and model4 from the tutorial (i.e., both for the cleaning & recode, as well as the linear model formulas). Do a modelsummary table displaying both model outputs. Interpret the coefficients for the MLR without the interaction (i.e., the first displayed model), and interpret model fit metrics for both models.

Loading Germany Data

germany_data <- read_fst("germany_data.fst")

Creating ‘weight’ column to account for weights for Germany

germany_data$weight <- germany_data$dweight * germany_data$pweight # creating weight column based on multiplication of d and p weight
survey_design <- svydesign(ids = ~1, data = germany_data, 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 what

Creating ‘auth’ variable

germany_data_model <- germany_data |> # now calculating authoritarian score
  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 )) # applying reverse coding - flips coding so that 6 = 1, 1 = 6 and adding flipped values up and scaling to 100

# Now you can calculate 'schwartzauth' after the NA recoding
germany_data_model$auth <- scales::rescale(germany_data_model$behave + 
                      germany_data_model$secure + 
                      germany_data_model$safety + 
                      germany_data_model$tradition + 
                      germany_data_model$rules, to=c(0,100), na.rm=TRUE)


germany_data_model <- germany_data_model |> filter(!is.na(auth)) # filtering NA

Double-checking if ‘auth’ variable was cleaned

unique(germany_data_model$auth)

##  [1]  16  84  88  56  60  92 100  96  76  68  44  40  72  80  64  36  48  24  52
## [20]  32  28  20   8  12   0   4

Re-coding ‘polID’ and creating ‘religious’ variable

germany_data_model <- germany_data_model |>
  mutate(
    polID = case_when( # creating new variable polID with lrscale
      lrscale %in% 0:3 ~ "Left",                    
      lrscale %in% 7:10 ~ "Right",                     
      lrscale %in% 4:6 ~ "Moderate",                  
      lrscale %in% c(77, 88, 99) ~ NA_character_      
    ),
   religious = case_when( # creating new variable = religious
      rlgdgr %in% c(77, 88, 99) ~ NA_real_,
      TRUE ~ rlgdgr
    )
  )

Double-checking if ‘polID’ and ‘religious’ variables were cleaned

unique(germany_data_model$polID)

## [1] "Right"    "Left"     "Moderate" NA

unique(germany_data_model$religious)

##  [1]  0  3  4  1  7  5  8  6 10  2  9 NA

germany_data_model <- germany_data_model |> filter(!is.na(polID)) # filtering NA
germany_data_model <- germany_data_model |> filter(!is.na(religious))

Creating models

model3 <- lm(auth ~ polID + religious, data = germany_data_model, weights = weight) # auth score, political ID and religious - created 2 models to try to account for variation in political score from previous models
model4 <- lm(auth ~ polID + religious + polID*religious, data = germany_data_model, weights = weight) # also creating interaction - is interaction between auth scores and political ID affected by religiosity

Creating modelsummary table

modelsummary(
  list(model3, model4),
  fmt = 1,
  estimate  = c( "{estimate} ({std.error}){stars}",
                "{estimate} ({std.error}){stars}"),
  statistic = NULL)

	(1)	(2)
(Intercept)	55.3 (0.3)***	56.4 (0.3)***
polIDModerate	4.6 (0.3)***	3.4 (0.4)***
polIDRight	7.7 (0.4)***	4.2 (0.7)***
religious	0.9 (0.0)***	0.6 (0.1)***
polIDModerate × religious		0.3 (0.1)***
polIDRight × religious		0.8 (0.1)***
Num.Obs.	23373	23373
R2	0.048	0.050
R2 Adj.	0.048	0.050
AIC	200594.1	200558.6
BIC	200634.4	200615.0
Log.Lik.	−100292.049	−100272.298
RMSE	17.27	17.24

Interpretation of Co-efficients

All the coefficients are significant with ***, indicating a p < 0.001. Meaning: There is less than a 0.1% probability that the observed result is due to chance and we incorrectly reject the null, under frequentist assumptions.

polIDModerate (4.6): When comparing individuals with the same religiousness, the average outcome for individuals identifying as politically moderate are 4.6 units higher compared to the reference category of politically left on the authoritarian values scale.
polIDRight (7.7): When comparing individuals with the same religiousness, individuals identifying as politically right are expected to be, on average, 7.7 units higher than the reference category of politically left on the authoritarian values scale.
religious (0.9): For the ‘religious’ predictor, a coefficient of 0.9 means that with every one-unit increase in religiousness, there is an associated expected average increase of 0.9 units on the authoritarian values scale, assuming that we are comparing people with the same political orientation.

Interpretation of Fit Metrics

Adjusted R-Squared (0.048 vs. 0.050): Model 1 explains 4.8% of variance in the outcome (auth scores) based on our religion and polID variables, which in comparison to model 2 that explains 5% is slightly less useful.

AIC (200594.1 vs. 200558.6): Model 1 has a higher AIC than model 2 with a 35.5 point difference, which means model 2 is a better descriptor of models with same dataset, as it has a lower AIC of 200558.6, indicating slightly better fit

BIC (200634.4 vs. 200615.0): Model 1 also has a higher BIC, which penalizes for any additional parameters added. Model 1 has a higher BIC than model 2 at 200615, suggesting model 1 has slightly more parameters it was penalized for and that model 2 is a slightly better fit.

Log.Lik. (-100292.049 vs. -100272.298): Model 1 has higher Log.Lik. at -100292.049 score than model 2 at -100272.298, indicating that model 2 is a better fit

RMSE (17.27 vs. 17.24): Finally, model one is 0.03 points higher than model 2, suggesting average magnitude of residuals is slightly higher for model 1, so model 2 is a better fit

Task 2

Now generate the model4 interaction plot that we did in the tutorial, but again using the German data instead of the French. Interpret.

Creating interaction plot

result <- predict_response(model4, terms = c("religious", "polID"))
plot(result)

Interpretation

There is a slight interaction effect of whether you are right leaning and religious to authoritarian score than to moderate and left, which suggests that there is a bit of a difference in religiosity and authoritarian scores if you identify as right leaning in Germany. As religiosity increases, authoritarian scores for politically oriented right individuals increases at a steeper slope, reaching almost 75 at maximum religiosity on our scale meaning that the lines are not parallel. Right oriented individuals remain the highest for authoritarian scores as religiosity increases. The interaction effect is generally low because ‘polID’ has a stronger effect on authoritarian scores, as seen by the coefficients on the table than if a person was moderate, right or left and religious.

Task 3

Use data for the Netherlands and the variable coding for model5 and model6 from the tutorial (i.e., both for the cleaning & recode, as well as the linear model formulas). Do a modelsummary table displaying both model outputs. Interpret the coefficients for the MLR without the interaction (i.e., the first displayed model), and interpret model fit metrics for both models.

Loading Netherlands Data

netherlands_data <- read_fst("netherlands_data.fst")

Creating ‘weight’ column to account for weights for Netherlands

netherlands_data$weight <- netherlands_data$dweight * netherlands_data$pweight # creating weight column based on multiplication of d and p weight
survey_design <- svydesign(ids = ~1, data = netherlands_data, 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 what

netherlands_data_model <- netherlands_data

Re-coding ‘gen’ variable

netherlands_data_model <- netherlands_data_model |>
  mutate(
    cohort = ifelse(yrbrn < 1930 | yrbrn > 2000, NA, yrbrn),
    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(yrbrn)  
    ),
    gen = factor(gen,
                 levels = c("1", "2", "3", "4"),
                 labels = c("Interwar", "Baby Boomers", "Gen X", "Millennials"))
  )

table(netherlands_data_model$gen)

## 
##     Interwar Baby Boomers        Gen X  Millennials 
##         3967         6359         4684         2803

table(netherlands_data_model$cohort)

## 
## 1930 1931 1932 1933 1934 1935 1936 1937 1938 1939 1940 1941 1942 1943 1944 1945 
##  125  155  143  140  154  177  190  160  201  218  229  216  231  275  280  253 
## 1946 1947 1948 1949 1950 1951 1952 1953 1954 1955 1956 1957 1958 1959 1960 1961 
##  378  356  332  298  325  285  325  292  394  311  375  336  295  336  324  328 
## 1962 1963 1964 1965 1966 1967 1968 1969 1970 1971 1972 1973 1974 1975 1976 1977 
##  358  361  350  310  342  362  364  385  373  358  345  282  273  250  273  284 
## 1978 1979 1980 1981 1982 1983 1984 1985 1986 1987 1988 1989 1990 1991 1992 1993 
##  250  233  219  198  218  203  204  202  213  192  179  150  142  143  122  117 
## 1994 1995 1996 1997 1998 1999 2000 
##  104  105   92   93   80   75   63

Coding Authoritarian Scale

netherlands_data_model <- netherlands_data_model |> # now calculating authoritarian score
  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 )) # applying reverse coding - flips coding so that 6 = 1, 1 = 6 and adding flipped values up and scaling to 100

# Now you can calculate 'schwartzauth' after the NA recoding
netherlands_data_model$auth <- scales::rescale(netherlands_data_model$behave + 
                      netherlands_data_model$secure + 
                      netherlands_data_model$safety + 
                      netherlands_data_model$tradition + 
                      netherlands_data_model$rules, to=c(0,100), na.rm=TRUE)


netherlands_data_model <- netherlands_data_model |> filter(!is.na(auth)) # filtering NA

unique(netherlands_data_model$auth)

##  [1]  80  72  48  36  56  64  52  28  76  60  68  92  44  84   8 100  40  88  32
## [20]  96  24  20  16   0   4  12

Re-coding ‘religion’

netherlands_data_model <- netherlands_data_model |>
  mutate(religion = case_when( # graph = different if change religious variable to categorical
    rlgblg == 2 ~ "No",
    rlgblg == 1 ~ "Yes",
    rlgblg %in% c(7, 8, 9) ~ NA_character_,
    TRUE ~ as.character(rlgblg)
  ))

# check
table(netherlands_data_model$religion)

## 
##    No   Yes 
## 10913  6750

Double-checking if ‘religion’ was cleaned

unique(netherlands_data_model$religion)

## [1] "No"  "Yes" NA

netherlands_data_model <- netherlands_data_model |> filter(!is.na(religion))

unique(netherlands_data_model$religion)

## [1] "No"  "Yes"

Re-coding ‘lrscale’

netherlands_data_model <- netherlands_data_model %>% # same with lrscale
  mutate(ID = case_when(
    lrscale >= 0 & lrscale <= 4 ~ "Left",
    lrscale >= 6 & lrscale <= 10 ~ "Right",
    lrscale > 10 ~ NA_character_,  # Set values above 10 as NA
    TRUE ~ NA_character_  # Ensure value 5 and any other unexpected values are set as NA
  ))
table(netherlands_data_model$ID)

## 
##  Left Right 
##  5439  7003

Double-checking if ‘ID’ was cleaned

unique(netherlands_data_model$ID)

## [1] "Left"  "Right" NA

netherlands_data_model <- netherlands_data_model |> filter(!is.na(ID))

unique(netherlands_data_model$ID)

## [1] "Left"  "Right"

Creating Linear Regressions

model5 <- lm(auth ~ religion + ID + gen, data = netherlands_data_model, weights = weight) # but with new religion and ID variables that were converted to categorical variables and including generation which is also categorical
model6 <- lm(auth ~ religion + ID + gen + religion*gen, data = netherlands_data_model, weights = weight) # same thing but interacting religion and generation

Creating modelsummary output

modelsummary(
  list(model5, model6),
  fmt = 1,
  estimate  = c( "{estimate} ({std.error}){stars}",
                "{estimate} ({std.error}){stars}"),
  statistic = NULL,
  coef_omit = "Intercept")

	(1)	(2)
religionYes	7.2 (0.3)***	5.7 (0.7)***
IDRight	4.2 (0.3)***	4.3 (0.3)***
genBaby Boomers	−5.2 (0.4)***	−5.9 (0.6)***
genGen X	−7.7 (0.4)***	−8.8 (0.6)***
genMillennials	−8.4 (0.5)***	−9.7 (0.6)***
religionYes × genBaby Boomers		1.1 (0.8)
religionYes × genGen X		2.3 (0.9)**
religionYes × genMillennials		2.9 (1.0)**
Num.Obs.	12116	12116
R2	0.124	0.125
R2 Adj.	0.124	0.125
AIC	101219.3	101213.5
BIC	101271.2	101287.5
Log.Lik.	−50602.672	−50596.757
RMSE	15.27	15.26

Interpretation of Co-efficients

religionYes (7.2): When comparing individuals with the same political orientation, the average authoritarian score for individuals who are religious are 7.2 units higher compared to the reference category of not religious.
IDRight (4.2): When comparing individuals with the same religiousness, individuals identifying as politically right are expected to be, on average, 4.2 units higher than the reference category of politically left on the authoritarian values scale.
genBaby Boomers (-5.2): When comparing individuals with the same religiousness, Baby Boomers are expected to be, on average, 5.2 units lower than the reference category of Interwar gen on the authoritarian values scale.
genGen X (-7.7): When comparing individuals with the same religiousness, Gen X are expected to be, on average, 7.7 units lower than the reference category of Interwar gen on the authoritarian values scale.
genMillennials (-8.4): When comparing individuals with the same religiousness, Millennials are expected to be, on average, 8.4 units lower than the reference category of Interwar gen on the authoritarian values scale. Millennials show the strongest decline.

Interpretation of Fit Metrics

Adjusted R-Squared (0.124 vs. 0.125): Model 1 explains 12.4% of variance in the outcome (auth scores) based on our religion and polID variables, which in comparison to model 2 that explains 12.5% is slightly less useful.

AIC (101219.3 vs. 101213.5): Model 1 has a higher AIC than model 2 with a 5.8 point difference, which means model 2 is a better descriptor of models with same dataset, as it has a lower AIC of 101213.5, indicating slightly better fit.

BIC (101271.2 vs. 101287.5): Model 1 also has a lower BIC at 101271.2, which penalizes for any additional parameters added. It suggests that model 1 has slightly good fit.

Log.Lik. (-50602.672 vs. -50596.757): Model 2 has lower Log.Lik. at -50596.757 than model 1, indicating that model 2 is a better fit.

RMSE (15.27 vs. 15.26): Finally, model 1 is 0.01 points higher than model 2, suggesting average magnitude of residuals is slightly higher for model 1, so model 2 is overall a better fit.

Task 4

Produce the model7 interaction plot from the tutorial, but again using the Netherlands data. Interpret.

netherlands_data_model_changed <- netherlands_data_model |>
  mutate(
    polID = case_when( # creating new variable polID with lrscale
      lrscale %in% 0:3 ~ "Left",                    
      lrscale %in% 7:10 ~ "Right",                     
      lrscale %in% 4:6 ~ "Moderate",                  
      lrscale %in% c(77, 88, 99) ~ NA_character_      
    ),
  )

unique(netherlands_data_model_changed$polID)

## [1] "Left"     "Right"    "Moderate"

Creating model7 plot

model7 <- lm(auth ~ + cohort + polID + cohort*polID, data = netherlands_data_model_changed, weights = weight)
interact_plot(model7, pred = cohort, modx = polID, jnplot = TRUE)

Interpretation

There is a stronger interaction effect between identifying as right and moderate for political orientation and authoritarian scores for each cohort in Netherlands because each line isn’t parallel. As each cohort progresses, authoritarian scores all decline, however, the decline is steeper for if you identify as right and moderate leaning than for those who identify on the left, with the right experiencing the sharpest decrease. It almost intersects with moderate scores, suggesting cohort has a strong effect on the positions of left and right leaning.

Task 5

Produce the model8 interaction plot from the tutorial, but again using the Netherlands data. Interpret.

Re-coding ‘religious’

netherlands_data_model_changed <- netherlands_data_model_changed |>
  mutate(
   religious = case_when( # creating new variable = religious
      rlgdgr %in% c(77, 88, 99) ~ NA_real_,
      TRUE ~ rlgdgr
    )
  )

Visualizing model8

model8 <- lm(auth ~ religious + ID + religious*ID, data = netherlands_data_model_changed, weights = weight)
interact_plot(model8, pred = religious, modx = ID, jnplot = TRUE)

Interpretation

There is also a slight interaction effect between identifying as left on the lrscale, religiosity and authoritarian scores in Netherlands. As you increase religiosity, those who identify both as politically left and right increase authoritarian scores, however, the scores of those who are right leaning are slightly steeper than those who identify on the left in Netherlands. Overall, the those who are left leaning have lower authoritarian scores than those who are right leaning, despite both experiencing an increase as the level of religiosity increases.

Nawaz_Fatima_Homework_5

Fatima Nawaz

2024-03-19

Loading Packages

Task 1

Loading Germany Data

Creating ‘weight’ column to account for weights for Germany

Creating ‘auth’ variable

Double-checking if ‘auth’ variable was cleaned

Re-coding ‘polID’ and creating ‘religious’ variable

Double-checking if ‘polID’ and ‘religious’ variables were cleaned

Creating models

Creating modelsummary table

Interpretation of Co-efficients

Interpretation of Fit Metrics

Task 2

Creating interaction plot

Interpretation

Task 3

Loading Netherlands Data

Creating ‘weight’ column to account for weights for Netherlands

Re-coding ‘gen’ variable

Coding Authoritarian Scale

Re-coding ‘religion’

Double-checking if ‘religion’ was cleaned

Re-coding ‘lrscale’

Double-checking if ‘ID’ was cleaned

Creating Linear Regressions

Creating modelsummary output

Interpretation of Co-efficients

Interpretation of Fit Metrics

Task 4

Creating model7 plot

Interpretation

Task 5

Re-coding ‘religious’

Visualizing model8

Interpretation