2. Main Text Figures
Figure 1: Conceptual Model and Data Flow
Figure 1: Conceptual framework of income-life satisfaction relationship across countries. This figure illustrates the theoretical pathways connecting individual income and life satisfaction, highlighting potential cultural and economic moderators that might influence this relationship across different national contexts.
Note: This figure is to be created externally as a vector graphic (180mm wide) and inserted into the manuscript. The conceptual framework integrates elements from economic theory (e.g., diminishing marginal utility), social comparison models, and cross-cultural psychology perspectives to explain how income-satisfaction associations might vary globally.
cat("Figure 1: Conceptual Model - To be created externally (180mm wide) and inserted into the manuscript.\n")## Figure 1: Conceptual Model - To be created externally (180mm wide) and inserted into the manuscript.Figure 2: Geographical Distribution of Average Life Satisfaction
Figure 2: Global distribution of average life satisfaction scores by country. This choropleth map displays mean life satisfaction scores (0-10 scale) for all countries included in the analysis, with darker blue indicating higher levels of reported life satisfaction.
Note: Data were aggregated from multiple longitudinal aging studies including CHARLS (China), ELSA (England), HRS (USA), KLOSA (Korea), MHAS (Mexico), and SHARE (Europe). Country-level averages were calculated using harmonized life satisfaction measures from individual responses collected between 2015 and 2023. Variations in average satisfaction demonstrate substantial cross-national differences that may be partially explained by cultural and economic factors examined in subsequent analyses.
# Calculate height in mm for saving function
fig2_height_mm_calc <- FIGURE_WIDTH_MM * FIG2_HEIGHT_FACTOR
# --- Data Preparation (same as before) ---
avg_satlife_fig2 <- clean_data %>% filter(!is.na(coun) & !is.na(satlife)) %>% group_by(coun) %>% summarise(mean_satlife = mean(satlife, na.rm = TRUE), .groups = 'drop') %>% filter(!is.na(mean_satlife))
world_map_sf_fig2_raw <- rnaturalearth::ne_countries(scale = "medium", returnclass = "sf") %>% select(iso_a3, name_en = name, sovereignt, geometry, area_proxy = scalerank) %>% filter(!is.na(iso_a3)) %>% sf::st_make_valid()
wintri_proj <- "+proj=wintri +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"
world_map_projected_fig2 <- world_map_sf_fig2_raw %>% sf::st_transform(crs = wintri_proj)
graticules_sf <- sf::st_graticule(ndiscr = 100, lat = seq(-90, 90, 15), lon = seq(-180, 180, 30)) %>% sf::st_transform(crs = wintri_proj)
world_map_joined_fig2 <- world_map_projected_fig2 %>% left_join(avg_satlife_fig2 %>% rename(mean_satlife_data = mean_satlife), by = c("iso_a3" = "coun"))
label_data_candidates_sf <- world_map_joined_fig2 %>% filter(!is.na(mean_satlife_data))
manual_label_iso <- c("CHN", "USA", "MEX", "KOR", "IND", "RUS", "JPN", "CAN", "GBR", "DEU", "FRA", "ESP", "ITA")
label_data_manual_sf <- label_data_candidates_sf %>% filter(iso_a3 %in% manual_label_iso)
label_data_repel_sf <- label_data_candidates_sf %>% filter(!iso_a3 %in% manual_label_iso, area_proxy <= 5)
label_data_manual_coords <- data.frame()
if(nrow(label_data_manual_sf) > 0) { centroids_manual <- sf::st_point_on_surface(label_data_manual_sf$geometry); temp_coords <- data.frame(iso_a3 = label_data_manual_sf$iso_a3, st_coordinates(centroids_manual)) %>% filter(!is.na(X) & !is.na(Y)); kor_offset_x <- 700000; kor_offset_y <- 150000; if ("KOR" %in% temp_coords$iso_a3) { temp_coords <- temp_coords %>% mutate(X = ifelse(iso_a3 == "KOR", X + kor_offset_x, X), Y = ifelse(iso_a3 == "KOR", Y + kor_offset_y, Y)) }; label_data_manual_coords <- temp_coords }
label_data_repel_coords <- data.frame()
if(nrow(label_data_repel_sf) > 0) { centroids_repel <- sf::st_point_on_surface(label_data_repel_sf$geometry); label_data_repel_coords <- data.frame(iso_a3 = label_data_repel_sf$iso_a3, st_coordinates(centroids_repel)) %>% filter(!is.na(X) & !is.na(Y)) }
min_lat_crop <- 0; max_lat_crop <- 85
full_xlim_proj <- st_bbox(world_map_projected_fig2)[c("xmin", "xmax")]
y_min_transformed <- st_sfc(st_point(c(0, min_lat_crop)), crs=4326) %>% st_transform(wintri_proj) %>% st_coordinates() %>% .[, "Y"]
y_max_transformed <- st_sfc(st_point(c(0, max_lat_crop)), crs=4326) %>% st_transform(wintri_proj) %>% st_coordinates() %>% .[, "Y"]
# --- Plot Creation (same as before, but using Blues palette) ---
plot_fig2 <- ggplot() +
geom_sf(data = graticules_sf, color = "gray85", linewidth = 0.15, linetype = "dotted") +
geom_sf(data = world_map_joined_fig2, aes(fill = mean_satlife_data), color = "gray50", linewidth = 0.12) +
{ if(nrow(label_data_manual_coords) > 0) { geom_text(data = label_data_manual_coords, aes(x = X, y = Y, label = iso_a3), size = 2.0, color = "black", fontface = "bold", check_overlap = FALSE) }} +
{ if(nrow(label_data_repel_coords) > 0) { ggrepel::geom_text_repel(data = label_data_repel_coords, aes(x = X, y = Y, label = iso_a3), size = 1.8, color = "black", fontface = "bold", min.segment.length = 0.1, box.padding = unit(0.15, "lines"), point.padding = unit(0.1, "lines"), segment.color = "gray30", segment.size = 0.20, force = 2.0, force_pull = 0.5, max.overlaps = Inf, bg.color = "white", bg.r = 0.08) }} +
scale_fill_distiller(palette = "Blues", direction = 1, name = "Mean Life Satisfaction (0-10)", limits = c(0, 10), na.value = "gray92", guide = guide_colorbar(barwidth = unit(FIGURE_WIDTH_MM * 0.30, "mm"), barheight = unit(3.0, "mm"), title.position = "top", title.hjust = 0.5, ticks.colour = "gray20", ticks.linewidth = 0.25, nbin=10)) +
coord_sf(crs = wintri_proj, xlim = full_xlim_proj, ylim = c(y_min_transformed, y_max_transformed), expand = FALSE, datum = NA) +
theme_void(base_family = "sans", base_size = 10) +
theme(legend.position = "bottom", legend.direction = "horizontal", legend.background = element_rect(fill="white", colour = NA), legend.title = element_text(size = rel(0.8), color = "black", face = "bold"), legend.text = element_text(size = rel(0.7), color = "black"), plot.background = element_rect(fill = "white", colour = NA), plot.margin = margin(t=1, r=1, b=1, l=1, unit="mm"))
# --- Save and Display ---
if (nrow(avg_satlife_fig2) > 0) {
# Save both EPS (for journal) and PNG (for HTML preview)
png_path <- save_journal_figure(plot_object = plot_fig2,
filename_base = "Fig2_Country_Life_Satisfaction",
height_mm = fig2_height_mm_calc,
width_mm = FIGURE_WIDTH_MM,
is_bitmap = FALSE) # 使用EPS格式，因为这是线图
# Display PNG in HTML output
knitr::include_graphics(png_path)
} else {
cat("Figure 2: No valid life satisfaction data found to display.\n")
plot_fig2_no_data <- ggplot() + theme_void() + annotate("text", x=0, y=0, label="Figure 2: No data available") + theme_publication(); print(plot_fig2_no_data)
}## Saved EPS: figures/Fig2_Country_Life_Satisfaction.eps (W: 180 mm, H: 99 mm)## Saved PNG: figures/Fig2_Country_Life_Satisfaction.png (W: 180 mm, H: 99 mm)
Figure 3: Meta-Analysis of Income-Life Satisfaction Associations
Figure 3: Forest plot of standardized income-life satisfaction effect sizes across countries. This meta-analysis displays standardized beta coefficients showing the relationship between income and life satisfaction for each analyzed country, along with 95% confidence intervals and an overall pooled estimate.
Note: Beta coefficients represent standardized effect sizes from regression models adjusting for age, gender, health status, and partnership status. GDP values are per-capita in international dollars (PPP adjusted). ATK represents the Atkinson inequality index (higher values indicate greater inequality). Sample size with analytical weight is shown for each country. The pooled estimate was derived using a random-effects model. Most countries show significant positive associations between income and life satisfaction, though with substantial cross-national heterogeneity in effect size magnitude.
# Check for required objects from the loaded RData file
if (!exists("model_meta_base") || !is.null(model_meta_base) && inherits(model_meta_base, "try-error")) {
cat("Error: 'model_meta_base' object for meta-analysis not found or is an error object. Forest plot cannot be generated.\n")
plot_fig3_placeholder <- ggplot() + theme_void() + annotate("text", x=0, y=0, label="Figure 3: Data missing (model_meta_base)")
print(plot_fig3_placeholder)
} else if (!exists("country_coefs_std_meta")) {
cat("Error: 'country_coefs_std_meta' (country-level coefficients for meta-analysis) not found. Forest plot cannot be generated.\n")
plot_fig3_placeholder <- ggplot() + theme_void() + annotate("text", x=0, y=0, label="Figure 3: Data missing (country_coefs_std_meta)")
print(plot_fig3_placeholder)
} else if (!exists("country_specific_models_output")) {
cat("Error: 'country_specific_models_output' (for model types) not found. Forest plot cannot be generated.\n")
plot_fig3_placeholder <- ggplot() + theme_void() + annotate("text", x=0, y=0, label="Figure 3: Data missing (country_specific_models_output)")
print(plot_fig3_placeholder)
} else {
# Prepare data for the forest plot
# Using 'country_coefs_std_meta' as the base, as it's loaded from RData
data_for_forest <- country_coefs_std_meta %>%
# Join model_type information
left_join(
country_specific_models_output %>% select(coun, model_type_from_source = model_type),
by = "coun"
) %>%
mutate(
# Use the 'coun' variable directly for study labels in the plot.
study_label = coun,
model_type = ifelse(!is.na(model_type_from_source), as.character(model_type_from_source), "N/A")
) # Removed arrange(estimate) to align with order="obs" from new snippet
# Calculate weights for display (as percentage)
# This uses 'w_perc_display' consistent with previous logic in satlife_visualizations_fix
data_for_forest$w_perc_display <- ""
w_percentage_values <- tryCatch(
weights(model_meta_base, type = "summary"),
error = function(e) {
warning("Could not retrieve weights from model_meta_base. N-based proportions will be used for display if N is available.")
return(NULL)
}
)
if (!is.null(w_percentage_values) && length(w_percentage_values) == nrow(data_for_forest)) {
data_for_forest$w_perc_display <- sprintf("%.1f", w_percentage_values)
} else if ("n_obs" %in% names(data_for_forest) && sum(data_for_forest$n_obs, na.rm = TRUE) > 0) {
total_n_obs <- sum(data_for_forest$n_obs, na.rm = TRUE)
if (total_n_obs > 0) {
data_for_forest$w_perc_display <- sprintf("%.1f", (data_for_forest$n_obs / total_n_obs) * 100)
}
}
# Use 'w_perc_display' for the text
n_weight_text_display <- paste0(
data_for_forest$n_obs,
ifelse(data_for_forest$w_perc_display != "" & !is.na(data_for_forest$w_perc_display),
paste0(" (", data_for_forest$w_perc_display, "%)"), # Using w_perc_display
"")
)
# Prepare ilab dataframe based on the new snippet's columns
ilab_vars <- data.frame(
GDP = ifelse(!is.na(data_for_forest$gdp), round(data_for_forest$gdp, 0), "N/A"),
ATK = ifelse(!is.na(data_for_forest$atk), formatC(data_for_forest$atk, digits = 2, format = "f"), "N/A"),
N_Weight = n_weight_text_display, # Using the calculated n_weight_text_display
stringsAsFactors = FALSE
)
# Set column name for N_Weight as per new snippet (without '%')
colnames(ilab_vars) <- c("GDP", "ATK", "N (Weight)")
# ilab positions from the new snippet
ilab_xpos <- c(-0.9, -0.7, -0.55)
# Dynamically calculate plot height
num_studies_fig3 <- nrow(data_for_forest)
# Adjusted spacing based on new snippet's implicit needs (fewer ilab, potentially different cex)
header_footer_space_mm <- 50
mm_per_study <- 6 # Adjusted from 5.5, can be tuned
fig3_height_mm <- header_footer_space_mm + (num_studies_fig3 * mm_per_study)
fig3_height_in <- mm_to_inches(fig3_height_mm)
fig3_width_in <- mm_to_inches(FIGURE_WIDTH_MM)
# Function to draw the forest plot, adapted from new snippet and existing structure
draw_forest_plot_fig3 <- function() {
par(family = "sans", mar=c(5,1,1.5,1)) # Margins from new snippet
summary_model <- summary(model_meta_base) # Keep for p-value check
summary_diamond_label <- "Overall Pooled Estimate" # Default label
if (model_meta_base$method == "FE") {
summary_diamond_label <- "Fixed-Effects Model"
} else if (model_meta_base$method %in% c("REML", "DL", "HE", "HS", "SJ", "PM", "EB")) {
summary_diamond_label <- "Random-Effects Model"
}
metafor::forest(model_meta_base, # Use the model from satlife_visualizations_fix
slab = data_for_forest$study_label, # 'coun' codes
ilab = ilab_vars,
ilab.xpos = ilab_xpos,
ilab.pos = 4,
ilab.lab = colnames(ilab_vars), # Use ilab.lab from new snippet
xlim = c(-1.4, 0.8), # From new snippet
alim = c(-0.1, 0.24), # From new snippet
refline = 0,
addpred = TRUE, # Maintained from satlife_visualizations_fix
predstyle = "bar", # From new snippet
pred.lab = summary_diamond_label, # Kept for clarity
showweights = FALSE,
digits = 3,
order = "obs", # From new snippet (data_for_forest not sorted by estimate)
psize = 1.5, # From new snippet
efac = c(1, 1), # From new snippet
cex = 0.8, # From new snippet
cex.lab = 1, # From new snippet
cex.axis = 0.85, # From new snippet
xlab = "Standardized Effect Size", # From new snippet
header = c("Country", "Effect [95% CI]"),
col = "darkblue",
border = "darkblue",
colout = "black", # From new snippet
shade = "zebra"
)
# Manual text for ilab headers is removed as ilab.lab is used.
# Corrected if condition for significance note (though mtext is commented out)
if (!is.null(summary_model$pval) &&
length(summary_model$pval) >= 1 && # Check if pval exists and has at least one element
!is.na(summary_model$pval[1]) && # Check if the first p-value is not NA
summary_model$pval[1] < 0.05) {
# Example mtext if re-enabled:
# mtext(paste0("Note: ", summary_diamond_label, " is statistically significant (p = ",
# format.pval(summary_model$pval[1], digits=3, eps=0.001), ")."),
# side=1, line=3.5, adj=0, cex=0.8*0.9) # Adjusted cex based on new base cex
}
}
# 对图3使用EPS格式（期刊要求的矢量格式）
eps_file_fig3 <- file.path("figures", "Fig3_Meta_Forest_Plot.eps")
postscript(file = eps_file_fig3, width = fig3_width_in, height = fig3_height_in,
horizontal = FALSE, paper = "special", family = "sans")
draw_forest_plot_fig3()
dev.off()
cat(paste("Saved forest plot as EPS:", eps_file_fig3, "\n"))
# 同时保存PNG用于HTML预览
png_file_fig3 <- file.path("figures", "Fig3_Meta_Forest_Plot.png")
png(png_file_fig3, width = FIGURE_WIDTH_MM, height = fig3_height_mm, units = "mm", res = 300, bg = "white")
draw_forest_plot_fig3()
dev.off()
cat(paste("Saved forest plot as PNG:", png_file_fig3, "\n"))
# Display PNG in HTML output
knitr::include_graphics(png_file_fig3)
}## Saved forest plot as EPS: figures/Fig3_Meta_Forest_Plot.eps## Saved forest plot as PNG: figures/Fig3_Meta_Forest_Plot.png
Forest plot of income-life satisfaction associations.
Figure 4: Correlation Matrix of Country-Level Variables
Figure 4: Correlation matrix of key country-level economic and cultural indicators. This matrix displays Pearson correlation coefficients between national-level variables used in moderator analyses, with scatter plots and fitted lines in the upper triangle and correlation values in the lower triangle.
Note: GNI and GDP represent per-capita values (PPP-adjusted international dollars). Atk (Atkinson index) and GINI are inequality measures. Cultural dimensions include: PDI (Power Distance), IDV (Individualism), MAS (Masculinity), UAI (Uncertainty Avoidance), LTO (Long-Term Orientation), and IVR (Indulgence vs. Restraint) from Hofstede’s cultural framework. Strong positive correlations exist between GDP and GNI (r = 0.95), and between IDV and GDP/GNI (r = 0.72/0.74), while strong negative correlations exist between inequality metrics and wealth indicators. Cultural dimensions show moderate intercorrelations, supporting their use as distinct moderators in subsequent analyses.
# Calculate height in mm for saving function
fig4_height_mm_calc <- FIGURE_WIDTH_MM * FIG4_HEIGHT_FACTOR
# --- Data Checks ---
if (!exists("country_level_desc_stats_data") || is.null(country_level_desc_stats_data) || nrow(country_level_desc_stats_data) == 0) {
cat("Data for Fig 4 (country_level_desc_stats_data) not found or empty.\n")
plot_fig4_placeholder <- ggplot() + theme_void() + annotate("text", x=0, y=0, label="Fig 4: Data missing") + theme_publication()
print(plot_fig4_placeholder)
} else {
# --- Data Preparation ---
fig4_filename_base <- "Fig4_Country_Level_Corr"
formal_var_names_fig4 <- c(gni = "GNI", gdp = "GDP", atk = "Atk", gini = "GINI", pdi = "PDI", idv = "IDV", mas = "MAS", uai = "UAI", lto = "LTO", ivr = "IVR")
data_for_corr_plot_renamed <- country_level_desc_stats_data %>% select(any_of(names(formal_var_names_fig4))) %>% rename_with(~ formal_var_names_fig4[.], .cols = any_of(names(formal_var_names_fig4)))
# --- Drawing Function ---
draw_describe_plot_func <- function() {
# No need to set par(family="sans") here, save_journal_figure handles it for base plots
bruceR::Describe(data_for_corr_plot_renamed, all.as.numeric = TRUE, digits = 2, plot = TRUE, upper.triangle = TRUE, upper.smooth = "lm")
}
# --- Save and Display ---
# 使用EPS格式保存（矢量图，期刊要求）
png_path <- save_journal_figure(plot_object = NULL,
filename_base = fig4_filename_base,
width_mm = FIGURE_WIDTH_MM,
height_mm = fig4_height_mm_calc,
is_ggplot = FALSE,
base_plot_draw_func = draw_describe_plot_func,
is_bitmap = FALSE) # 使用EPS格式
# Display PNG in HTML output
if (!is.null(png_path)) {
knitr::include_graphics(png_path)
} else {
cat("Warning: PNG file for Fig 4 could not be created.\n")
}
}## Descriptive Statistics:
## ──────────────────────────────────────────────────────────────────────────
## N Mean SD | Median Min Max Skewness Kurtosis
## ──────────────────────────────────────────────────────────────────────────
## GNI 33 36392.12 13540.29 | 34920.00 13140.00 69430.00 0.58 -0.30
## GDP 33 37862.91 17267.38 | 35220.98 13201.67 107888.58 1.97 5.76
## Atk 33 0.18 0.06 | 0.17 0.11 0.36 0.91 0.39
## GINI 33 0.33 0.05 | 0.33 0.25 0.48 0.72 0.66
## PDI 33 51.85 21.28 | 50.00 11.00 100.00 0.14 -0.55
## IDV 33 56.42 20.22 | 60.00 18.00 91.00 -0.37 -0.99
## MAS 33 47.73 23.00 | 47.00 5.00 100.00 0.05 -0.63
## UAI 33 70.76 21.90 | 74.00 23.00 112.00 -0.44 -0.52
## LTO 33 57.64 18.76 | 58.00 24.00 100.00 0.16 -0.70
## IVR 33 46.45 21.90 | 48.00 13.00 97.00 0.20 -1.06
## ──────────────────────────────────────────────────────────────────────────## Saved EPS: figures/Fig4_Country_Level_Corr.eps (W: 180 mm, H: 153 mm)
## Descriptive Statistics:
## ──────────────────────────────────────────────────────────────────────────
## N Mean SD | Median Min Max Skewness Kurtosis
## ──────────────────────────────────────────────────────────────────────────
## GNI 33 36392.12 13540.29 | 34920.00 13140.00 69430.00 0.58 -0.30
## GDP 33 37862.91 17267.38 | 35220.98 13201.67 107888.58 1.97 5.76
## Atk 33 0.18 0.06 | 0.17 0.11 0.36 0.91 0.39
## GINI 33 0.33 0.05 | 0.33 0.25 0.48 0.72 0.66
## PDI 33 51.85 21.28 | 50.00 11.00 100.00 0.14 -0.55
## IDV 33 56.42 20.22 | 60.00 18.00 91.00 -0.37 -0.99
## MAS 33 47.73 23.00 | 47.00 5.00 100.00 0.05 -0.63
## UAI 33 70.76 21.90 | 74.00 23.00 112.00 -0.44 -0.52
## LTO 33 57.64 18.76 | 58.00 24.00 100.00 0.16 -0.70
## IVR 33 46.45 21.90 | 48.00 13.00 97.00 0.20 -1.06
## ──────────────────────────────────────────────────────────────────────────## Saved PNG: figures/Fig4_Country_Level_Corr.png (W: 180 mm, H: 153 mm)
Figure 5: Moderation by Economic and Inequality Indicators
Figure 5: Economic and inequality indicators as moderators of income-life satisfaction associations. These scatter plots show how country-level economic factors modify the relationship between income and life satisfaction, with each point representing a country’s standardized effect size and the size of points proportional to precision (1/SE).
Note: All moderators are standardized (z-scores). Panels show: (A) GDP per capita, (B) GNI per capita, (C) Atkinson inequality index, and (D) GINI coefficient. Regression lines show fitted meta-regression models weighted by effect size precision. Countries with higher GDP and GNI per capita show significantly stronger positive associations between income and life satisfaction (p < 0.01). Inequality indicators (Atkinson index and GINI) demonstrate a negative moderation effect, with higher inequality associated with weaker income-satisfaction relationships, suggesting that the psychological benefits of income may be diminished in more unequal societies.
# Calculate height in mm for saving function
fig5_height_mm_calc <- FIGURE_WIDTH_MM * FIG5_HEIGHT_FACTOR
# --- Data Checks ---
if (!exists("country_coefs_std_meta")) {
cat("Data for Fig 5 (country_coefs_std_meta) not found.\n")
plot_A <- ggplot() + theme_void() + labs(title="A") + theme_publication()
plot_B <- ggplot() + theme_void() + labs(title="B") + theme_publication()
plot_C <- ggplot() + theme_void() + labs(title="C") + theme_publication()
plot_D <- ggplot() + theme_void() + labs(title="D") + theme_publication()
plot_fig5_combined <- (plot_A | plot_B) / (plot_C | plot_D)
print(plot_fig5_combined)
} else {
# --- Data Preparation & Plotting ---
moderators_fig5_def <- list("gdp_z" = "A: GDP per Capita (z)", "gni_z" = "B: GNI per Capita (z)", "atk_z" = "C: Atkinson Index (z)", "gini_z" = "D: GINI Coefficient (z)")
plot_list_fig5 <- list()
moderator_data_fig5 <- country_coefs_std_meta %>% mutate(precision_weight = 1 / sqrt(vi))
for (mod_var_fig5_idx in seq_along(moderators_fig5_def)) {
mod_var_name_fig5 <- names(moderators_fig5_def)[mod_var_fig5_idx]
full_mod_title <- moderators_fig5_def[[mod_var_name_fig5]]
subplot_letter_fig5 <- stringr::str_extract(full_mod_title, "^[A-D]")
x_axis_title_fig5 <- stringr::str_replace(full_mod_title, "^[A-D]:\\s+", "")
if (mod_var_name_fig5 %in% names(moderator_data_fig5)) {
num_distinct_mod_fig5 <- moderator_data_fig5 %>% filter(!is.na(.data[[mod_var_name_fig5]])) %>% pull(.data[[mod_var_name_fig5]]) %>% n_distinct()
p_fig5_panel <- ggplot(moderator_data_fig5, aes_string(x = mod_var_name_fig5, y = "estimate")) +
geom_point(aes(size = precision_weight), alpha = 0.6, color = palette_discrete_viridis[1]) +
ggrepel::geom_text_repel(aes(label = coun), size = 2.2, max.overlaps = Inf, segment.alpha = 0.5, segment.color = "grey40", min.segment.length = 0, box.padding = 0.25, point.padding = 0.2, bg.color = "white", bg.r = 0.05) +
scale_size_continuous(guide = "none") +
labs(title = subplot_letter_fig5, x = x_axis_title_fig5, y = if(mod_var_fig5_idx %% 2 == 1) expression("Country Effect Size ("*beta*")") else "") +
theme_publication()
if (num_distinct_mod_fig5 >= 2) {
p_fig5_panel <- p_fig5_panel + geom_smooth(method = "lm", aes(weight = precision_weight), color = "black", se = TRUE, linewidth=0.7)
}
plot_list_fig5[[mod_var_name_fig5]] <- p_fig5_panel
} else {
plot_list_fig5[[mod_var_name_fig5]] <- ggplot() + annotate("text",x=0,y=0,label=paste(mod_var_name_fig5, "\nnot found.")) + theme_void() + labs(title=subplot_letter_fig5, x=x_axis_title_fig5) + theme_publication()
}
}
# --- Combine and Save ---
if (length(plot_list_fig5) == 4) {
plot_fig5_combined <- wrap_plots(plot_list_fig5, ncol = 2)
# 使用EPS格式保存（矢量图，期刊要求）
png_path <- save_journal_figure(plot_fig5_combined,
"Fig5_Eco_Ineq_Moderators",
height_mm = fig5_height_mm_calc,
width_mm = FIGURE_WIDTH_MM,
is_bitmap = FALSE) # 使用EPS格式
# Display PNG in HTML output
if (!is.null(png_path)) {
knitr::include_graphics(png_path)
} else {
cat("Warning: PNG file for Fig 5 could not be created.\n")
}
} else {
cat("Could not generate all 4 panels for Figure 5.\n")
plot_fig5_combined_error <- ggplot() + theme_void() + labs(title="Figure 5 Error") + theme_publication()
print(plot_fig5_combined_error)
}
}## Saved EPS: figures/Fig5_Eco_Ineq_Moderators.eps (W: 180 mm, H: 171 mm)## Saved PNG: figures/Fig5_Eco_Ineq_Moderators.png (W: 180 mm, H: 171 mm)
Figure 6: Moderation by Cultural Dimensions
Figure 6: Hofstede’s cultural dimensions as moderators of income-life satisfaction associations. These scatter plots illustrate how each of Hofstede’s six cultural dimensions modifies the relationship between income and life satisfaction, with each point representing a country’s standardized effect size and point size proportional to precision (1/SE).
Note: All cultural dimensions are standardized (z-scores). Panels show: (A) Power Distance (PDI), (B) Individualism (IDV), (C) Masculinity (MAS), (D) Uncertainty Avoidance (UAI), (E) Long-Term Orientation (LTO), and (F) Indulgence vs. Restraint (IVR). Regression lines show fitted meta-regression models weighted by effect size precision. Most notably, individualism (IDV) shows a significant positive moderation effect (p < 0.01), with more individualistic societies exhibiting stronger relationships between income and life satisfaction. Power distance (PDI) shows a significant negative moderation effect (p < 0.05), suggesting that hierarchical societies may experience decreased psychological benefits from personal income. Other cultural dimensions show non-significant moderation trends that warrant further investigation with larger country samples.
# Calculate height in mm for saving function
fig6_height_mm_calc <- FIGURE_WIDTH_MM * FIG6_HEIGHT_FACTOR
# --- Data Checks ---
if (!exists("country_coefs_std_meta")) {
cat("Data for Fig 6 (country_coefs_std_meta) not found.\n")
plot_list_fig6_init <- replicate(6, ggplot() + theme_void() + labs(title="Data Missing") + theme_publication(), simplify = FALSE)
plot_fig6_combined <- wrap_plots(plot_list_fig6_init, ncol=2)
print(plot_fig6_combined)
} else {
# --- Data Preparation & Plotting ---
moderators_fig6_def <- list("pdi_z" = "A: Power Distance (PDI, z)", "idv_z" = "B: Individualism (IDV, z)", "mas_z" = "C: Masculinity (MAS, z)", "uai_z" = "D: Uncertainty Avoidance (UAI, z)", "lto_z" = "E: Long-Term Orientation (LTO, z)", "ivr_z" = "F: Indulgence vs. Restraint (IVR, z)")
plot_list_fig6 <- list()
moderator_data_fig6 <- country_coefs_std_meta %>% mutate(precision_weight = 1 / sqrt(vi))
for (i in seq_along(moderators_fig6_def)) {
mod_var_name_fig6 <- names(moderators_fig6_def)[i]
full_mod_title_fig6 <- moderators_fig6_def[[mod_var_name_fig6]]
subplot_letter_fig6 <- stringr::str_extract(full_mod_title_fig6, "^[A-F]")
x_axis_title_fig6 <- stringr::str_replace(full_mod_title_fig6, "^[A-F]:\\s+", "")
if (mod_var_name_fig6 %in% names(moderator_data_fig6)) {
num_distinct_mod_fig6 <- moderator_data_fig6 %>% filter(!is.na(.data[[mod_var_name_fig6]])) %>% pull(.data[[mod_var_name_fig6]]) %>% n_distinct()
p_fig6_panel <- ggplot(moderator_data_fig6, aes_string(x = mod_var_name_fig6, y = "estimate")) +
geom_point(aes(size = precision_weight), alpha = 0.6, color = viridis::viridis_pal(option="plasma")(8)[i+1]) +
ggrepel::geom_text_repel(aes(label = coun), size = 2.2, max.overlaps = Inf, segment.alpha = 0.5, segment.color = "grey50", min.segment.length = 0, box.padding = 0.25, point.padding = 0.2, bg.color = "white", bg.r = 0.05) +
scale_size_continuous(guide = "none") +
labs(title = subplot_letter_fig6, x = x_axis_title_fig6, y = if(i %% 2 == 1) expression("Country Effect Size ("*beta*")") else "") +
theme_publication()
if (num_distinct_mod_fig6 >= 2) {
p_fig6_panel <- p_fig6_panel + geom_smooth(method = "lm", aes(weight = precision_weight), color = "black", se = TRUE, linewidth=0.7)
}
plot_list_fig6[[mod_var_name_fig6]] <- p_fig6_panel
} else {
plot_list_fig6[[mod_var_name_fig6]] <- ggplot() + annotate("text",x=0,y=0,label=paste(mod_var_name_fig6, "\nnot found.")) + theme_void() + labs(title=subplot_letter_fig6, x=x_axis_title_fig6) + theme_publication()
}
}
# --- Combine and Save ---
if (length(plot_list_fig6) == 6) {
plot_fig6_combined <- wrap_plots(plot_list_fig6, ncol = 2, nrow = 3)
# 使用EPS格式保存（矢量图，期刊要求）
png_path <- save_journal_figure(plot_fig6_combined,
"Fig6_Cultural_Moderators",
height_mm = fig6_height_mm_calc,
width_mm = FIGURE_WIDTH_MM,
is_bitmap = FALSE) # 使用EPS格式
# Display PNG in HTML output
if (!is.null(png_path)) {
knitr::include_graphics(png_path)
} else {
cat("Warning: PNG file for Fig 6 could not be created.\n")
}
} else {
cat("Could not generate all 6 panels for Figure 6.\n")
plot_fig6_combined_error <- ggplot() + theme_void() + labs(title="Figure 6 Error") + theme_publication()
print(plot_fig6_combined_error)
}
}## Saved EPS: figures/Fig6_Cultural_Moderators.eps (W: 180 mm, H: 243 mm)## Saved PNG: figures/Fig6_Cultural_Moderators.png (W: 180 mm, H: 243 mm)