As a reminder, the study 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
install.packages("tidyverse", repos = "https://cloud.r-project.org/")
The downloaded binary packages are in
/var/folders/4v/p0llgt5n3_l5x55x1g2swx_h0000gn/T//RtmpmRwQ7x/downloaded_packageslibrary(tidyverse)
Group7 <- read_csv("Desktop/Group7.csv")
Group7 <- Group7 %>%
filter(ideo != 8)#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)#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#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")
The downloaded binary packages are in
/var/folders/4v/p0llgt5n3_l5x55x1g2swx_h0000gn/T//RtmpmRwQ7x/downloaded_packageslibrary(modelsummary)
modelsummary(list("Regression Results" = regressionmodel),
output = "kableExtra", # Can also use "gt", "html", "latex"
stars = TRUE, statistic = "p.value",
title = "Linear Regression Results")| Regression Results | |
|---|---|
| (Intercept) | 4.068*** |
| (<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. Our model predicts that for every one unit increase in age, we predict a 0.110 increase in units of happiness, and that for every one unit increase in educational level, we predict a 0.146 unit increase in units of 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 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!).
#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 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")))
# Table of happiness values by religiosity levels and political orientation
library(dplyr)
library(knitr)
library(kableExtra)
Group7 %>%
filter(
!is.na(political_orientation_categorical),
!is.na(religiosity_categorical),
!is.na(happy)
) %>%
group_by(political_orientation_categorical, religiosity_categorical) %>%
summarise(
n = n(),
mean_happy = mean(happy)
) %>%
ungroup() %>%
kable(
digits = 2,
caption = "Mean Happiness and Sample Size by Political Orientation and Religiosity"
) %>%
kable_styling(full_width = FALSE, bootstrap_options = c("striped", "hover"))| political_orientation_categorical | religiosity_categorical | n | mean_happy |
|---|---|---|---|
| Conservative | High | 35 | 5.26 |
| Conservative | Low | 21 | 4.62 |
| Conservative | Medium | 25 | 4.88 |
| Liberal | High | 27 | 5.78 |
| Liberal | Low | 25 | 4.60 |
| Liberal | Medium | 26 | 5.00 |
| Moderate | High | 39 | 5.23 |
| Moderate | Low | 54 | 4.39 |
| Moderate | Medium | 54 | 4.52 |
# Plot code
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 = "Happiness", linetype = "Political Orientation", title = "Mean Happiness Levels by Religiosity and Political Orientation") +
theme_minimal() +
theme(legend.position = c(0.25, 0.85), # Adjust coordinates to fine-tune placement
legend.direction = "vertical",
legend.background = element_rect(fill = "white", color = NA), # Optional: add a background for clarity
legend.text = element_text(size = 9),
legend.title = element_text(size = 11)) +
scale_linetype_manual(values = c("dashed", "dotted", "solid"),
labels = c("Conservative", "Liberal", "Moderate"))
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).