1. Average Level of Concern for Each Societal Threat

This analysis examines the average level of concern for different societal threats across all participants.

# Load required libraries
library(ggplot2)
library(dplyr)
library(tidyr)
library(readr)

# Read in the dataset - adjust the path to match your file location
data <- read.csv("/Users/linjinhu/Desktop/Jackson Research Assistant Data Task 2/Dataset.csv")

# Create a new data frame for societal threats
societal_threats <- data %>%
  select(starts_with("SocietalThreat")) %>%
  gather(key = "Threat", value = "ConcernLevel") %>%
  mutate(
    Threat = recode(Threat,
                    "SocietalThreat_1" = "Natural Disaster",
                    "SocietalThreat_2" = "Famine",
                    "SocietalThreat_3" = "Pollution",
                    "SocietalThreat_4" = "Attack by a terrorist group",
                    "SocietalThreat_5" = "Crime surge",
                    "SocietalThreat_6" = "Discrimination",
                    "SocietalThreat_7" = "Increase in National Debt",
                    "SocietalThreat_8" = "Surge of COVID-19 cases",
                    "SocietalThreat_9" = "Influx of refugees",
                    "SocietalThreat_10" = "Illegal immigration")
  )

# Calculate mean and standard error for each threat
summary_stats <- societal_threats %>%
  group_by(Threat) %>%
  summarise(
    Mean = mean(ConcernLevel, na.rm = TRUE),
    SE = sd(ConcernLevel, na.rm = TRUE) / sqrt(n())
  )

# Create the plot with error bars
ggplot(summary_stats, aes(x = reorder(Threat, Mean), y = Mean)) +
  geom_bar(stat = "identity", fill = "skyblue", color = "black", width = 0.7) +
  geom_errorbar(aes(ymin = Mean - SE, ymax = Mean + SE), width = 0.2) +
  coord_flip() +  # Rotate the plot for better readability
  labs(
    title = "Average Level of Concern for Each Societal Threat",
    x = "Societal Threat",
    y = "Average Concern Level",
    caption = "Note: Error bars represent the standard error of the mean"
  ) +
  theme_minimal() +
  theme(
    text = element_text(size = 12),
    plot.title = element_text(hjust = 0.5, size = 16, face = "bold"),
    axis.title = element_text(size = 14),
    axis.text = element_text(size = 12)
  )

2. Reliability Analysis of Cultural Tightness Scale

# Load necessary libraries
library(psych)  

# Specify the cultural tightness variables
ct_variables <- c("CT_1", "CT_2", "CT_3", "CT_4", "CT_5", "CT_6", "CT_7", "CT_8", "CT_9")

# Calculate Cronbach's alpha for the cultural tightness scale
ct_reliability <- psych::alpha(data[ct_variables])

# Display the reliability results
print(ct_reliability$total$raw_alpha)
## [1] 0.852059
print(ct_reliability$total)
##  raw_alpha std.alpha   G6(smc) average_r     S/N         ase     mean       sd
##   0.852059 0.8518593 0.8433304 0.3898446 5.75034 0.002365539 4.463244 1.114378
##   median_r
##  0.3825554
# Check item-total correlations and alpha if item deleted
print(ct_reliability$alpha.drop)
##      raw_alpha std.alpha   G6(smc) average_r      S/N    alpha se       var.r
## CT_1 0.8455301 0.8453422 0.8335414 0.4059063 5.465888 0.002489429 0.005468475
## CT_2 0.8465528 0.8467426 0.8357902 0.4085015 5.524971 0.002479941 0.005261906
## CT_3 0.8449542 0.8444038 0.8338900 0.4041808 5.426893 0.002497249 0.006043600
## CT_4 0.8368918 0.8363429 0.8243047 0.3897944 5.110335 0.002622785 0.006784724
## CT_5 0.8299128 0.8300920 0.8151420 0.3791490 4.885539 0.002750163 0.004598951
## CT_6 0.8300108 0.8302190 0.8162206 0.3793610 4.889940 0.002743616 0.004747957
## CT_7 0.8281025 0.8276720 0.8160154 0.3751410 4.802888 0.002771282 0.005420764
## CT_8 0.8359187 0.8355867 0.8245712 0.3884835 5.082232 0.002643620 0.006495159
## CT_9 0.8296732 0.8294524 0.8180409 0.3780837 4.863467 0.002745450 0.006260900
##          med.r
## CT_1 0.3930975
## CT_2 0.3979028
## CT_3 0.3943025
## CT_4 0.3780364
## CT_5 0.3825554
## CT_6 0.3780364
## CT_7 0.3743173
## CT_8 0.3780364
## CT_9 0.3617490
# Get overall descriptive statistics
mean_score <- mean(rowMeans(data[ct_variables], na.rm = TRUE), na.rm = TRUE)
sd_score <- sd(rowMeans(data[ct_variables], na.rm = TRUE), na.rm = TRUE)

print(paste("Mean score:", round(mean_score, 2)))
## [1] "Mean score: 4.46"
print(paste("Standard deviation:", round(sd_score, 2)))
## [1] "Standard deviation: 1.11"

The 9-item Cultural Tightness scale demonstrated good internal consistency reliability (Cronbach’s α = .85). The mean score across all participants was 4.46 (SD = 1.11) on a 7-point scale. Item analysis revealed that all items contributed positively to the scale’s reliability, with the alpha-if-item-deleted values ranging from .83 to .85, indicating that no single item’s removal would substantially improve the scale’s overall reliability. Item-total correlations were adequate across all items. These results suggest that the desire for cultural tightness scale is a reliable measure of the construct in this sample, exceeding the conventional threshold of .70 for acceptable reliability in social science research.

3. Correlations Between Cultural Tightness and Societal Threats

# Load necessary libraries
library(dplyr)
library(knitr)
library(kableExtra)

# Define Fisher's r-to-z transformation function for calculating CIs
fisher.r.to.z <- function(r, n) {
  z <- 0.5 * log((1 + r) / (1 - r))
  se <- 1 / sqrt(n - 3)
  z.lower <- z - 1.96 * se
  z.upper <- z + 1.96 * se
  r.lower <- (exp(2 * z.lower) - 1) / (exp(2 * z.lower) + 1)
  r.upper <- (exp(2 * z.upper) - 1) / (exp(2 * z.upper) + 1)
  return(list(lower = r.lower, upper = r.upper))
}

# Calculate mean Cultural Tightness score for each participant
ct_variables <- c("CT_1", "CT_2", "CT_3", "CT_4", "CT_5", "CT_6", "CT_7", "CT_8", "CT_9")
data$CT_mean <- rowMeans(data[ct_variables], na.rm = TRUE)

# Define threat variables and their labels
threat_vars <- paste0("SocietalThreat_", 1:10)
threat_labels <- c("Natural Disaster", "Famine", "Pollution", 
                   "Attack by a Terrorist Group", "Crime Surge", 
                   "Discrimination", "Increase in National Debt", 
                   "A Surge of COVID-19 Cases", "Influx of Refugees", 
                   "Illegal Immigration")

# Create a data frame to store correlation results
corr_table <- data.frame(
  Threat = threat_labels,
  r = numeric(10),
  p_value = numeric(10),
  n = numeric(10),
  ci_lower = numeric(10),
  ci_upper = numeric(10),
  stringsAsFactors = FALSE
)

# Calculate correlations between CT_mean and each threat variable
for (i in 1:10) {
  # Count valid pairs (non-missing values)
  valid_n <- sum(!is.na(data$CT_mean) & !is.na(data[[threat_vars[i]]]))
  
  # Calculate correlation
  corr_result <- cor.test(data$CT_mean, data[[threat_vars[i]]], 
                          use = "pairwise.complete.obs")
  
  corr_table$r[i] <- corr_result$estimate
  corr_table$p_value[i] <- corr_result$p.value
  corr_table$n[i] <- valid_n
  
  # Calculate 95% CI using Fisher's r-to-z transformation
  ci <- fisher.r.to.z(corr_table$r[i], valid_n)
  corr_table$ci_lower[i] <- ci$lower
  corr_table$ci_upper[i] <- ci$upper
}

# Round correlation coefficients and CIs to 2 decimal places (APA standard)
corr_table$r <- round(corr_table$r, 2)
corr_table$ci_lower <- round(corr_table$ci_lower, 2)
corr_table$ci_upper <- round(corr_table$ci_upper, 2)

# Add significance stars
corr_table$sig <- ""
corr_table$sig[corr_table$p_value < 0.05] <- "*"
corr_table$sig[corr_table$p_value < 0.01] <- "**"
corr_table$sig[corr_table$p_value < 0.001] <- "***"

# Combine r with significance stars
corr_table$r_formatted <- paste0(corr_table$r, corr_table$sig)

# Create 95% CI column
corr_table$ci_formatted <- paste0("[", corr_table$ci_lower, ", ", corr_table$ci_upper, "]")

# Create final table for display
final_table <- corr_table[, c("Threat", "r_formatted", "ci_formatted")]
colnames(final_table) <- c("Societal Threat", "r", "95% CI")

# Render the table with kableExtra for enhanced formatting
final_table %>%
  kable(
    caption = "Correlations Between Desire for Cultural Tightness and Level of Concern for Societal Threats",
    format = "html",
    align = c("l", "c", "c"),
    escape = FALSE
  ) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = FALSE,
    position = "center"
  ) %>%
  add_header_above(c(" " = 1, "Correlation with Cultural Tightness" = 2)) %>%
  footnote(
    general = paste0(
      "N = ", min(corr_table$n), ". ",
      "Correlations represent the relationship between participants' desire for cultural tightness ",
      "(measured as the average of nine cultural tightness items) and their level of concern for each societal threat. ",
      "Higher positive correlations indicate that as desire for cultural tightness increases, ",
      "concern for the corresponding threat also tends to increase. ",
      "95% CI = 95% confidence intervals.",
      "\n* p < .05. ** p < .01. *** p < .001."
    ),
    general_title = "Note. ",
    footnote_order = "general"
  )
Correlations Between Desire for Cultural Tightness and Level of Concern for Societal Threats
Correlation with Cultural Tightness
Societal Threat r 95% CI
Natural Disaster 0.09*** [0.07, 0.11]
Famine 0.1*** [0.08, 0.12]
Pollution 0.11*** [0.09, 0.13]
Attack by a Terrorist Group 0.05*** [0.03, 0.07]
Crime Surge 0.2*** [0.18, 0.22]
Discrimination 0.08*** [0.06, 0.1]
Increase in National Debt 0.16*** [0.14, 0.18]
A Surge of COVID-19 Cases 0.11*** [0.09, 0.13]
Influx of Refugees 0.18*** [0.16, 0.2]
Illegal Immigration 0.17*** [0.15, 0.19]
Note.
N = 8688. Correlations represent the relationship between participants’ desire for cultural tightness (measured as the average of nine cultural tightness items) and their level of concern for each societal threat. Higher positive correlations indicate that as desire for cultural tightness increases, concern for the corresponding threat also tends to increase. 95% CI = 95% confidence intervals.
* p < .05. ** p < .01. *** p < .001.

4.Regression Analysis of COVID-19 Concern and Cultural Tightness

# Load necessary libraries
library(dplyr)
library(broom) 

# Calculate mean Cultural Tightness score for each participant if not already done
if(!"CT_mean" %in% names(data)) {
  ct_variables <- c("CT_1", "CT_2", "CT_3", "CT_4", "CT_5", "CT_6", "CT_7", "CT_8", "CT_9")
  data$CT_mean <- rowMeans(data[ct_variables], na.rm = TRUE)
}

# Identify the COVID-19 variable
covid_var <- "SocietalThreat_8" # A surge of COVID-19 cases

# Create a clean dataset for analysis
analysis_data <- data %>%
  select(CT_mean, !!sym(covid_var)) %>%
  rename(covid_concern = !!sym(covid_var)) %>%
  filter(!is.na(CT_mean), !is.na(covid_concern))

# Run regression model
covid_model <- lm(CT_mean ~ covid_concern, data = analysis_data)

# Get model summary
model_summary <- summary(covid_model)

# Get standardized beta coefficient (effect size)
sd_x <- sd(analysis_data$covid_concern)
sd_y <- sd(analysis_data$CT_mean)
std_beta <- model_summary$coefficients[2,1] * (sd_x / sd_y)

# Calculate Cohen's f² effect size
r_squared <- model_summary$r.squared
f_squared <- r_squared / (1 - r_squared)

# Calculate correlation as well (alternative effect size measure)
correlation <- cor(analysis_data$CT_mean, analysis_data$covid_concern, 
                   use = "pairwise.complete.obs")

# Get tidy version of model results
tidy_results <- tidy(covid_model, conf.int = TRUE)

# Show diagnostic information
print(model_summary)
## 
## Call:
## lm(formula = CT_mean ~ covid_concern, data = analysis_data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.5573 -0.5694 -0.0382  0.6650  2.8385 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    4.06255    0.04186  97.058   <2e-16 ***
## covid_concern  0.09894    0.00991   9.984   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.108 on 8686 degrees of freedom
## Multiple R-squared:  0.01135,    Adjusted R-squared:  0.01123 
## F-statistic: 99.68 on 1 and 8686 DF,  p-value: < 2.2e-16
cat("\n")
cat("Sample size:", nrow(analysis_data), "\n")
## Sample size: 8688
cat("Standardized beta coefficient:", round(std_beta, 3), "\n")
## Standardized beta coefficient: 0.107
cat("Cohen's f²:", round(f_squared, 3), "\n")
## Cohen's f²: 0.011
cat("Correlation:", round(correlation, 3), "\n")
## Correlation: 0.107
cat("95% CI for unstandardized coefficient:", 
    paste0("[", round(tidy_results$conf.low[2], 3), ", ", 
           round(tidy_results$conf.high[2], 3), "]"), "\n\n")
## 95% CI for unstandardized coefficient: [0.08, 0.118]

A linear regression analysis was conducted to examine the association between concern about a rise in COVID-19 cases and desire for cultural tightness. COVID-19 concern significantly predicted desire for cultural tightness, b = 0.099, 95% CI [0.08, 0.118], t(8686) = 9.98, p < .001. The standardized coefficient (β = 0.107) indicates that a one standard deviation increase in COVID-19 concern is associated with a 0.107 standard deviation increase in desire for cultural tightness. The model explains 1.1% of the variance in cultural tightness (R² = 0.011), representing a small effect size (Cohen’s f² = 0.011). This finding suggests that individuals who express greater concern about COVID-19 surges also tend to desire greater cultural tightness, though the magnitude of this association is relatively small.

5.Association Between COVID-19 Concern and Cultural Tightness Across Countries

library(ggplot2)
library(dplyr)

# Calculate mean Cultural Tightness score for each participant if not already done
if(!"CT_mean" %in% names(data)) {
  ct_variables <- c("CT_1", "CT_2", "CT_3", "CT_4", "CT_5", "CT_6", "CT_7", "CT_8", "CT_9")
  data$CT_mean <- rowMeans(data[ct_variables], na.rm = TRUE)
}

# Identify the COVID-19 variable
covid_var <- "SocietalThreat_8" # A surge of COVID-19 cases

# Select three countries for comparison
selected_countries <- c("Canada", "Columbia", "Japan")

# Filter data for the selected countries
plot_data <- data %>%
  filter(country %in% selected_countries) %>%
  filter(!is.na(CT_mean), !is.na(!!sym(covid_var))) %>%
  select(country, covid_concern = !!sym(covid_var), cultural_tightness = CT_mean)

# Calculate correlations for each country
correlations <- plot_data %>%
  group_by(country) %>%
  summarize(
    correlation = cor(covid_concern, cultural_tightness, use = "pairwise.complete.obs"),
    n = n(),
    .groups = "drop"
  ) %>%
  mutate(
    # Calculate p-values
    p_value = mapply(function(r, n) {
      t_value <- r * sqrt((n - 2) / (1 - r^2))
      2 * pt(-abs(t_value), df = n - 2)  # Two-tailed test
    }, correlation, n),
    # Add significance indicators
    significance = case_when(
      p_value < 0.001 ~ "***",
      p_value < 0.01 ~ "**",
      p_value < 0.05 ~ "*",
      TRUE ~ ""
    ),
    # Create labels
    label = sprintf("r = %.2f%s", correlation, significance)
  )

# Create a more readable country name for Columbia
plot_data <- plot_data %>%
  mutate(country_label = case_when(
    country == "Columbia" ~ "Colombia",
    TRUE ~ country
  ))

# Update correlations table with the same label
correlations <- correlations %>%
  mutate(country_label = case_when(
    country == "Columbia" ~ "Colombia",
    TRUE ~ country
  ))

# Create a clean, focused scatterplot
ggplot(plot_data, aes(x = covid_concern, y = cultural_tightness, color = country_label)) +
  # Add points with slight jitter to avoid overplotting
  geom_jitter(alpha = 0.6, width = 0.2, height = 0, size = 1.5) +
  # Add regression lines for each country
  geom_smooth(method = "lm", se = TRUE, aes(fill = country_label), alpha = 0.2) +
  # Add correlation annotations
  geom_text(data = correlations, 
            aes(x = 1.5, y = 6.8 - match(country, selected_countries) * 0.3, 
                label = label, color = country_label),
            hjust = 0, size = 4, fontface = "bold", show.legend = FALSE) +
  # Set proper axis scales
  scale_x_continuous(breaks = 1:5, limits = c(0.5, 5.5)) +
  scale_y_continuous(limits = c(1, 7)) +
  # Labels
  labs(
    title = "Association Between COVID-19 Concern and Cultural Tightness",
    x = "Concern About COVID-19 Surge (1 = Low, 5 = High)",
    y = "Desire for Cultural Tightness (1-7 Scale)",
    color = "Country",
    fill = "Country"
  ) +
  # Clean theme
  theme_light() +
  theme(
    legend.position = "bottom",
    panel.grid.minor = element_blank(),
    axis.title = element_text(face = "bold"),
    axis.text = element_text(size = 10),
    legend.text = element_text(size = 10)
  )

Interpretation of Cross-Country Comparison

The scatterplot reveals interesting variations in the relationship between COVID-19 concern and cultural tightness across the three countries. Japan shows little to no correlation, suggesting that in Japanese society, concern about COVID-19 is not meaningfully associated with desires for cultural tightness. Colombia shows a modest positive correlation, while Canada shows the most interesting pattern with a slight negative correlation—suggesting that in Canada, greater COVID-19 concern might actually be associated with slightly lower desire for cultural tightness. These varying patterns highlight the importance of considering cultural context when examining the relationship between societal threats and preferences for cultural norms.