PSY225D Replication Study Data Analysis

As a reminder, the studying we are aiming to replicate is this:

Bixter, M. T. (2015). Happiness, political orientation, and religiosity. Personality and Individual Differences, 72, 7-11. ://www.sciencedirect.com/science/article/pii/S0191886914004553

Code chunk for getting things set up/a tiny bit of cleaning on our political views variable:

Code 
install.packages("tidyverse", repos = "https://cloud.r-project.org/")
Installing package into '/Users/claireheyne/Library/R/arm64/4.3/library'
(as 'lib' is unspecified)

The downloaded binary packages are in
    /var/folders/4v/p0llgt5n3_l5x55x1g2swx_h0000gn/T//RtmpVrS60R/downloaded_packages
Code 
library(tidyverse)
── 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.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
Code 
Group7 <- read_csv("Desktop/Group7.csv")
Rows: 996 Columns: 28
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
dbl (28): gender, age, ethnic_1, ethnic_2, ethnic_3, ethnic_4, ethnic_5, eth...

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Code 
Group7 <- Group7 %>% 
  filter(ideo != 8)

Code chunk for standardizing our religiosity variables and creating our religiosity index:

Code 
#standardizing our religiosity variables
Group7 <- 
  Group7 %>% 
mutate(pray = as.numeric(g7_pray),  
         pray_s = scale(g7_pray))    

Group7 <- 
  Group7 %>% 
  mutate(attend = as.numeric(g7_attend),  
         attend_s = scale(g7_attend))    

Group7 <- 
  Group7 %>% 
  mutate(relition = as.numeric(g7_relition),  
         relition_s = scale(g7_relition))    

Group7 <- 
  Group7 %>% 
  mutate(god = as.numeric(g7_god),  
         god_s = scale(g7_god))    

Group7 <- 
  Group7 %>% 
  mutate(sprtprsn = as.numeric(g7_sprtprsn),  
         sprtprsn_s = scale(g7_sprtprsn))    

Group7 <- 
  Group7 %>% 
  mutate(relprsn = as.numeric(g7_relprsn),  
         relprsn_s = scale(g7_relprsn))    

Group7 <- 
  Group7 %>% 
  mutate(relaffil = as.numeric(relaffil),  
         relaffil_s = scale(relaffil))    

#creating the religiosity index

Group7 <-
  Group7 %>%
  mutate(religiosity = (pray_s + attend_s + relition_s + god_s + sprtprsn_s + relprsn_s + relaffil_s)/7)

Code chunk for our regression model, in which we predicted happiness levels based on political orientation, religiosity, age, sex, educational level, and the interaction between political orientation and religiosity.

Code 
#writing the code for our regression model

#first, we ensure r views all of our variables as numeric for the sake of the regression output

Group7$ideo <- as.numeric(as.character(Group7$ideo))  

Group7$religiosity <- as.numeric(as.character(Group7$religiosity))

Group7$age <- as.numeric(Group7$age)  

Group7$educ <- as.numeric(as.character(Group7$educ))

#regression code

regressionmodel = lm(formula = happy ~ ideo + religiosity + age + gender + educ + ideo:religiosity, data=Group7)

summary(regressionmodel)

Call:
lm(formula = happy ~ ideo + religiosity + age + gender + educ + 
    ideo:religiosity, data = Group7)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.1579 -0.9683  0.3935  1.1730  2.6549 

Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)       4.06826    0.43039   9.452  < 2e-16 ***
ideo             -0.05849    0.05150  -1.136  0.25696    
religiosity       0.61809    0.31555   1.959  0.05107 .  
age               0.11031    0.05573   1.979  0.04869 *  
gender           -0.01539    0.15841  -0.097  0.92265    
educ              0.14587    0.05113   2.853  0.00464 ** 
ideo:religiosity -0.01232    0.07290  -0.169  0.86587    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.578 on 297 degrees of freedom
  (608 observations deleted due to missingness)
Multiple R-squared:  0.09275,   Adjusted R-squared:  0.07442 
F-statistic:  5.06 on 6 and 297 DF,  p-value: 5.823e-05
Code 
#some code for making our regression results more well-formatted (showing intercepts and p-values in parentheses)

options(repos = c(CRAN = "https://cloud.r-project.org/"))
install.packages("modelsummary")
Installing package into '/Users/claireheyne/Library/R/arm64/4.3/library'
(as 'lib' is unspecified)

The downloaded binary packages are in
    /var/folders/4v/p0llgt5n3_l5x55x1g2swx_h0000gn/T//RtmpVrS60R/downloaded_packages
Code 
library(modelsummary)
Warning: package 'modelsummary' was built under R version 4.3.3
Code 
modelsummary(list("Regression Results" = regressionmodel),
             output = "kableExtra",  # Can also use "gt", "html", "latex"
             stars = TRUE, statistic = "p.value",
             title = "Linear Regression Results")
Linear Regression Results
Regression Results
(Intercept) 4.068***
(&lt;0.001)
ideo −0.058
(0.257)
religiosity 0.618+
(0.051)
age 0.110*
(0.049)
gender −0.015
(0.923)
educ 0.146**
(0.005)
ideo × religiosity −0.012
(0.866)
Num.Obs. 304
R2 0.093
R2 Adj. 0.074
AIC 1148.8
BIC 1178.6
Log.Lik. −566.410
RMSE 1.56
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001

So, to summarize, we found that age and education were significant predictors of happiness, such that greater age predicts higher levels of happiness and higher educational levels predict higher happiness.

Although we might say that our religiosity predictor was approaching statistical significance, it still was not a significant predictor of happiness, which differs from Bixter’s finding that higher levels of religiosity significantly predicted higher levels of happiness. Also contrary to what Bixter found, neither political orientation nor the interaction between political orientation and religiosity significantly predicted happiness in our model. In Bixter’s study, being more conservative was associated with greater happiness, and it was found that religiosity has a greater effect on happiness for more conservative individuals compared to those who are more politically liberal.

Our adjusted r-squared value of 0.074 indicates that this only 7.4% of the variance in happiness levels is explained by the variance in political orientation, religiosity, age, sex, educational level, and the interaction between political orientation and religiosity. This is actually slightly better than the R-squared values obtained by Bixter (0.04 and 0.02!).

Code chunk for a data visualization of mean happiness of different political orientations by religiosity level

Code 
#data visualization for mean happiness of different political orientations by religiosity level

#first, let's make political orientation a categorical variable

Group7 <- Group7 %>%
  mutate(political_orientation_categorical = case_when(
    ideo %in% c(1, 2) ~ "Liberal",
    ideo %in% c(3, 4, 5) ~ "Moderate",
    ideo %in% c(6, 7) ~ "Conservative",
    TRUE ~ NA_character_  
  ))

#and we will do the same for religiosity

Group7$religiosity <- as.numeric(as.character(Group7$religiosity))

#let's determine the 33rd and 67th percentiles so we can split the data into low, moderate, and high religiosity

  quantile(Group7$religiosity, probs = c(0.33, 0.67), na.rm = TRUE)
       33%        67% 
-0.2320644  0.3445929
Code 
Group7 <- Group7 %>%
  mutate(religiosity_categorical = case_when(
    religiosity < -0.2320644 ~ "Low",
    religiosity >= -0.2320644 & religiosity <= 0.3445929 ~ "Medium",
    religiosity > 0.3445929 ~ "High",
    TRUE ~ NA_character_  
  ))

Group7$happy <- as.numeric(as.character(Group7$happy))

data_summary <- Group7 %>%
  group_by(religiosity_categorical, political_orientation_categorical) %>%
  summarize(mean_happy = mean(happy, na.rm = TRUE), .groups = "drop")

data_summary <- data_summary %>%
  mutate(religiosity_categorical = factor(religiosity_categorical, levels = c("Low", "Medium", "High")))

# Plot

data_summary <- na.omit(data_summary)

data_summary %>%
ggplot(aes(x = religiosity_categorical, y = mean_happy,
                         group = political_orientation_categorical,
                         linetype = political_orientation_categorical)) +
  geom_point() +
  geom_line() +
  labs(x = "Religiosity", y = "Mean Happiness", linetype = "Political Orientation") +
  theme_minimal() +
  theme(legend.position = "top") +
  theme(legend.text = element_text(size = 9),     
      legend.title = element_text(size = 11)) + 
  scale_linetype_manual(values = c("solid", "longdash", "dotted"))

So, to summarize what we can see in the visual (this does NOT mean these observations are statistically significant!), being more religious is associated with greater levels of happiness, regardless of political orientation. Generally, liberals tend to be happier than conservatives, who tend to be happier than moderates. Highly religious liberals are the happiest group, although the low religiosity conservatives and liberals demonstrate a very similar level of mean happiness (and the high religiosity conservatives and moderates also demonstrate a very similar level of mean happiness).