1 OSF Repository Layout
2 Load Data
3 Variable Construction
4 Descriptives
5 Primary Confirmatory Models
6 Continuous-Age Simple Slopes
7 Interaction Figures
8 Bonferroni-Corrected Correlations
9 Robustness Checks
10 Session Information

1 OSF Repository Layout

Expected repository structure:

data/
  analysis_dataset.csv       # preferred de-identified analysis dataset
  OnlyGamers.xlsx            # optional raw/working dataset, if ethics permit
analysis/
  ACP_25_0797_OSF_Reproducible_Analysis.Rmd
outputs/
  tables/
  figures/

This document first looks for data/analysis_dataset.csv. If that file is absent, it attempts to read data/OnlyGamers.xlsx and rebuild the analysis variables from the GSQ items.

2 Load Data

analysis_csv <- "analysis_dataset.csv"
raw_xlsx <- "RawAdultGamers"

if (file.exists(analysis_csv)) {
  dat <- readr::read_csv(analysis_csv, show_col_types = FALSE)
  source_file <- analysis_csv
} else if (file.exists(raw_xlsx)) {
  dat <- readxl::read_excel(raw_xlsx)
  source_file <- raw_xlsx
} else {
  stop("No dataset found. Add data/analysis_dataset.csv or data/OnlyGamers.xlsx.")
}

cat("Source dataset:", source_file, "\n")

## Source dataset: analysis_dataset.csv

cat("Rows:", nrow(dat), "Columns:", ncol(dat), "\n")

## Rows: 507 Columns: 165

3 Variable Construction

zscore <- function(x) as.numeric(scale(as.numeric(x)))

make_age_group <- function(age) {
  dplyr::case_when(
    age <= 35 ~ "18-35",
    age <= 55 ~ "36-55",
    TRUE ~ "56+"
  )
}

genre_map <- tibble::tribble(
  ~code,  ~label,                    ~freq_item, ~ability_item,
  "SPG",  "Sports games",            "GSQ1",     "GSQ2",
  "FPSG", "First-person shooter",    "GSQ3",     "GSQ4",
  "RPG",  "Role-playing",            "GSQ5",     "GSQ6",
  "AAG",  "Action-adventure",        "GSQ7",     "GSQ8",
  "STG",  "Strategy",                "GSQ9",     "GSQ10",
  "PSG",  "Puzzle-solving",          "GSQ11",    "GSQ12"
)

dat <- dat %>%
  mutate(
    Age = as.numeric(Age),
    AgeGroup = if ("AgeGroup" %in% names(.)) AgeGroup else make_age_group(Age),
    AgeGroup = factor(AgeGroup, levels = c("18-35", "36-55", "56+")),
    Sex_Male = if ("Sex_Male" %in% names(.)) {
      as.numeric(Sex_Male)
    } else {
      as.numeric(str_to_lower(as.character(Sex)) == "male")
    },
    Education = as.numeric(Education)
  )

for (i in seq_len(nrow(genre_map))) {
  code <- genre_map$code[i]
  freq_item <- genre_map$freq_item[i]
  ability_item <- genre_map$ability_item[i]

  dat[[paste0(code, "_Freq")]] <- as.numeric(dat[[freq_item]])
  dat[[paste0(code, "_Ability")]] <- as.numeric(dat[[ability_item]])
  dat[[paste0(code, "_CurrentSkill")]] <-
    dat[[paste0(code, "_Freq")]] + dat[[paste0(code, "_Ability")]]
}

freq_cols <- paste0(genre_map$code, "_Freq")
ability_cols <- paste0(genre_map$code, "_Ability")
legacy_skill_cols <- paste0(genre_map$code, "_CurrentSkill")

dat <- dat %>%
  mutate(
    GSQ_Freq_Total_Rebuilt = rowSums(across(all_of(freq_cols)), na.rm = TRUE),
    GSQ_Ability_Total_Rebuilt = rowSums(across(all_of(ability_cols)), na.rm = TRUE),
    GSQ_Total_Skill_Rebuilt = rowSums(across(all_of(legacy_skill_cols)), na.rm = TRUE)
  )

# Legacy frequency group retained only for supplemental/continuity analyses.
# This reproduces the submitted manuscript's stratified median split by age group.
dat <- dat %>%
  group_by(AgeGroup) %>%
  mutate(
    Legacy_Frequency_Group = if_else(
      GSQ_Freq_Total_Rebuilt >= median(GSQ_Freq_Total_Rebuilt, na.rm = TRUE),
      "High",
      "Low"
    )
  ) %>%
  ungroup() %>%
  mutate(Legacy_Frequency_Group = factor(Legacy_Frequency_Group, levels = c("Low", "High")))

rebuild_checks <- tibble(
  legacy_GSQ_Frequency_matches_rebuilt =
    if ("GSQ_Frequency" %in% names(dat)) all(dat$GSQ_Frequency == dat$GSQ_Freq_Total_Rebuilt, na.rm = TRUE) else NA,
  legacy_GSQ_Ability_matches_rebuilt =
    if ("GSQ_Ability" %in% names(dat)) all(dat$GSQ_Ability == dat$GSQ_Ability_Total_Rebuilt, na.rm = TRUE) else NA,
  legacy_GSQ_Total_matches_rebuilt =
    if ("GSQ_Total" %in% names(dat)) all(dat$GSQ_Total == dat$GSQ_Total_Skill_Rebuilt, na.rm = TRUE) else NA
)

readr::write_csv(rebuild_checks, "outputs/tables/rebuild_checks.csv")
kable(rebuild_checks)

legacy_GSQ_Frequency_matches_rebuilt	legacy_GSQ_Ability_matches_rebuilt	legacy_GSQ_Total_matches_rebuilt
TRUE	TRUE	TRUE

continuous_to_z <- c(
  "Age", "Education",
  "GSQ_Freq_Total_Rebuilt", "GSQ_Ability_Total_Rebuilt",
  "GSQ_Total_Skill_Rebuilt",
  freq_cols, ability_cols, legacy_skill_cols
)

for (v in continuous_to_z) {
  dat[[paste0(v, "_z")]] <- zscore(dat[[v]])
}

dat <- dat %>%
  mutate(
    Age_z_x_GSQ_Freq_Total_z = Age_z * GSQ_Freq_Total_Rebuilt_z,
    Age_z_x_GSQ_Ability_Total_z = Age_z * GSQ_Ability_Total_Rebuilt_z,
    GSQ_Freq_Total_z_sq = GSQ_Freq_Total_Rebuilt_z^2,
    Age_z_x_GSQ_Freq_Total_z_sq = Age_z * GSQ_Freq_Total_z_sq,
    Age_z_x_GSQ_Total_Skill_z = Age_z * GSQ_Total_Skill_Rebuilt_z
  )

for (code in genre_map$code) {
  dat[[paste0("Age_z_x_", code, "_Freq_z")]] <-
    dat$Age_z * dat[[paste0(code, "_Freq_z")]]
  dat[[paste0("Age_z_x_", code, "_Ability_z")]] <-
    dat$Age_z * dat[[paste0(code, "_Ability_z")]]
}

readr::write_csv(dat, "outputs/tables/analysis_dataset_rebuilt.csv")

4 Descriptives

age_descriptives <- dat %>%
  count(AgeGroup, Sex, name = "n") %>%
  group_by(AgeGroup) %>%
  mutate(percent_within_age_group = 100 * n / sum(n)) %>%
  ungroup()

legacy_frequency_counts <- dat %>%
  count(AgeGroup, Legacy_Frequency_Group, name = "n") %>%
  group_by(AgeGroup) %>%
  mutate(percent_within_age_group = 100 * n / sum(n)) %>%
  ungroup()

primary_descriptives <- dat %>%
  summarise(
    n = n(),
    age_mean = mean(Age, na.rm = TRUE),
    age_sd = sd(Age, na.rm = TRUE),
    freq_mean = mean(GSQ_Freq_Total_Rebuilt, na.rm = TRUE),
    freq_sd = sd(GSQ_Freq_Total_Rebuilt, na.rm = TRUE),
    ability_mean = mean(GSQ_Ability_Total_Rebuilt, na.rm = TRUE),
    ability_sd = sd(GSQ_Ability_Total_Rebuilt, na.rm = TRUE),
    esq_mean = mean(Total_Score_ESQ, na.rm = TRUE),
    esq_sd = sd(Total_Score_ESQ, na.rm = TRUE),
    emq_mean = mean(`EMQ-R`, na.rm = TRUE),
    emq_sd = sd(`EMQ-R`, na.rm = TRUE),
    prmq_mean = mean(PRMQ, na.rm = TRUE),
    prmq_sd = sd(PRMQ, na.rm = TRUE)
  )

readr::write_csv(age_descriptives, "outputs/tables/descriptives_agegroup.csv")
readr::write_csv(legacy_frequency_counts, "outputs/tables/descriptives_legacy_frequency_group.csv")
readr::write_csv(primary_descriptives, "outputs/tables/descriptives_primary.csv")

kable(primary_descriptives, digits = 2)

n	age_mean	age_sd	freq_mean	freq_sd	ability_mean	ability_sd	esq_mean	esq_sd	emq_mean	emq_sd	prmq_mean	prmq_sd
507	42.74	13.89	13.02	5.66	17.89	5.74	20.44	11.4	23.72	9	35.49	11.05

kable(legacy_frequency_counts, digits = 2)

AgeGroup	Legacy_Frequency_Group	n	percent_within_age_group
18-35	Low	83	44.62
18-35	High	103	55.38
36-55	Low	95	43.58
36-55	High	123	56.42
56+	Low	42	40.78
56+	High	61	59.22

5 Primary Confirmatory Models

Primary outcomes:

ESQ-R total: Total_Score_ESQ
EMQ-R total: EMQ-R
PRMQ total: PRMQ

Primary predictors:

continuous age, z-scored
total self-reported GSQ frequency, z-scored
total perceived GSQ ability, z-scored
Age x Frequency
Age x Ability
sex and education covariates

HC3 robust standard errors are used. Holm correction is applied across the 15 focal confirmatory tests.

primary_outcomes <- c(
  "Total_Score_ESQ" = "ESQ-R total",
  "EMQ-R" = "EMQ-R total",
  "PRMQ" = "PRMQ total"
)

term_labels <- c(
  "Age_z" = "Age",
  "GSQ_Freq_Total_Rebuilt_z" = "GSQ self-reported frequency",
  "GSQ_Ability_Total_Rebuilt_z" = "GSQ perceived gaming ability",
  "Age_z_x_GSQ_Freq_Total_z" = "Age x Frequency",
  "Age_z_x_GSQ_Ability_Total_z" = "Age x Ability",
  "Sex_Male" = "Sex: male",
  "Education_z" = "Education"
)

primary_terms <- c(
  "Sex_Male",
  "Education_z",
  "Age_z",
  "GSQ_Freq_Total_Rebuilt_z",
  "GSQ_Ability_Total_Rebuilt_z",
  "Age_z_x_GSQ_Freq_Total_z",
  "Age_z_x_GSQ_Ability_Total_z"
)

focal_terms <- c(
  "Age_z",
  "GSQ_Freq_Total_Rebuilt_z",
  "GSQ_Ability_Total_Rebuilt_z",
  "Age_z_x_GSQ_Freq_Total_z",
  "Age_z_x_GSQ_Ability_Total_z"
)

fit_lm_hc3 <- function(data, outcome, terms) {
  formula <- as.formula(paste0("`", outcome, "` ~ ", paste(terms, collapse = " + ")))
  model <- lm(formula, data = data)
  vc <- sandwich::vcovHC(model, type = "HC3")
  ct <- lmtest::coeftest(model, vcov. = vc)
  ci <- coef(model) %>%
    names() %>%
    map_dfr(function(term) {
      est <- coef(model)[[term]]
      se <- sqrt(diag(vc))[term]
      tibble(
        term = term,
        estimate = est,
        robust_se = se,
        conf_low = est - 1.96 * se,
        conf_high = est + 1.96 * se
      )
    })
  broom::tidy(ct) %>%
    rename(
      b = estimate,
      statistic = statistic,
      p_raw = p.value
    ) %>%
    left_join(ci, by = "term") %>%
    mutate(
      b = estimate,
      r2 = summary(model)$r.squared,
      adj_r2 = summary(model)$adj.r.squared,
      n = stats::nobs(model),
      aic = AIC(model),
      bic = BIC(model)
    )
}

primary_model_results <- purrr::imap_dfr(primary_outcomes, function(label, outcome) {
  fit_lm_hc3(dat, outcome, primary_terms) %>%
    filter(term != "(Intercept)") %>%
    mutate(
      outcome = outcome,
      outcome_label = label,
      term_label = recode(term, !!!term_labels),
      family = if_else(term %in% focal_terms, "primary_focal", "covariate")
    )
})

primary_model_results <- primary_model_results %>%
  group_by(family) %>%
  mutate(
    p_holm = if_else(family == "primary_focal", p.adjust(p_raw, method = "holm"), NA_real_),
    p_fdr = if_else(family == "primary_focal", p.adjust(p_raw, method = "BH"), NA_real_)
  ) %>%
  ungroup()

readr::write_csv(primary_model_results, "outputs/tables/PrimaryModels_Coefficients_R.csv")

primary_table <- primary_model_results %>%
  transmute(
    Outcome = outcome_label,
    Predictor = term_label,
    b = round(b, 3),
    SE_HC3 = round(robust_se, 3),
    CI_low = round(conf_low, 3),
    CI_high = round(conf_high, 3),
    p_raw = signif(p_raw, 3),
    p_holm_focal = signif(p_holm, 3)
  )

kable(primary_table)

Outcome	Predictor	b	SE_HC3	CI_low	CI_high	p_raw	p_holm_focal
ESQ-R total	Sex: male	0.117	0.935	-1.716	1.949	9.01e-01	NA
ESQ-R total	Education	-1.896	0.441	-2.761	-1.032	2.06e-05	NA
ESQ-R total	Age	-4.756	0.491	-5.718	-3.794	0.00e+00	0.00e+00
ESQ-R total	GSQ self-reported frequency	3.338	0.620	2.124	4.552	1.00e-07	1.10e-06
ESQ-R total	GSQ perceived gaming ability	-3.227	0.553	-4.312	-2.143	0.00e+00	1.00e-07
ESQ-R total	Age x Frequency	3.363	0.681	2.028	4.698	1.10e-06	9.70e-06
ESQ-R total	Age x Ability	-1.695	0.603	-2.878	-0.512	5.17e-03	2.58e-02
EMQ-R total	Sex: male	0.259	0.837	-1.382	1.899	7.57e-01	NA
EMQ-R total	Education	-0.507	0.398	-1.287	0.274	2.04e-01	NA
EMQ-R total	Age	-3.194	0.401	-3.979	-2.408	0.00e+00	0.00e+00
EMQ-R total	GSQ self-reported frequency	2.956	0.521	1.936	3.976	0.00e+00	3.00e-07
EMQ-R total	GSQ perceived gaming ability	-1.694	0.426	-2.528	-0.860	7.94e-05	6.35e-04
EMQ-R total	Age x Frequency	1.554	0.619	0.340	2.768	1.24e-02	4.96e-02
EMQ-R total	Age x Ability	-0.322	0.452	-1.208	0.563	4.76e-01	4.76e-01
PRMQ total	Sex: male	-1.299	1.121	-3.496	0.898	2.47e-01	NA
PRMQ total	Education	-0.747	0.454	-1.637	0.142	1.00e-01	NA
PRMQ total	Age	-3.449	0.498	-4.425	-2.473	0.00e+00	0.00e+00
PRMQ total	GSQ self-reported frequency	2.982	0.775	1.463	4.500	1.35e-04	9.43e-04
PRMQ total	GSQ perceived gaming ability	-1.807	0.907	-3.585	-0.029	4.70e-02	1.39e-01
PRMQ total	Age x Frequency	2.553	0.715	1.152	3.955	3.92e-04	2.35e-03
PRMQ total	Age x Ability	-1.595	0.798	-3.159	-0.030	4.63e-02	1.39e-01

6 Continuous-Age Simple Slopes

linear_combo <- function(model, vcov_mat, weights) {
  beta <- coef(model)
  common <- intersect(names(weights), names(beta))
  w <- rep(0, length(beta))
  names(w) <- names(beta)
  w[common] <- weights[common]
  est <- sum(w * beta)
  se <- sqrt(as.numeric(t(w) %*% vcov_mat %*% w))
  z <- est / se
  p <- 2 * pnorm(abs(z), lower.tail = FALSE)
  tibble(b = est, SE = se, z = z, p = p)
}

probe_values <- c("-1 SD" = -1, "Mean" = 0, "+1 SD" = 1)

simple_slopes_agez <- purrr::imap_dfr(primary_outcomes, function(label, outcome) {
  formula <- as.formula(paste0("`", outcome, "` ~ ", paste(primary_terms, collapse = " + ")))
  model <- lm(formula, data = dat)
  vc <- sandwich::vcovHC(model, type = "HC3")

  purrr::imap_dfr(probe_values, function(age_z, probe_label) {
    freq <- linear_combo(
      model, vc,
      c(
        "GSQ_Freq_Total_Rebuilt_z" = 1,
        "Age_z_x_GSQ_Freq_Total_z" = age_z
      )
    ) %>%
      mutate(Predictor = "Frequency")

    ability <- linear_combo(
      model, vc,
      c(
        "GSQ_Ability_Total_Rebuilt_z" = 1,
        "Age_z_x_GSQ_Ability_Total_z" = age_z
      )
    ) %>%
      mutate(Predictor = "Ability")

    bind_rows(freq, ability) %>%
      mutate(
        Outcome = label,
        Age_probe = probe_label,
        Age_z = age_z
      )
  })
})

readr::write_csv(simple_slopes_agez, "outputs/tables/SimpleSlopes_AgeZ_R.csv")
kable(simple_slopes_agez, digits = 3)

b	SE	z	p	Predictor	Outcome	Age_probe	Age_z
-0.025	0.747	-0.034	0.973	Frequency	ESQ-R total	-1 SD	-1
-1.532	0.761	-2.015	0.044	Ability	ESQ-R total	-1 SD	-1
3.338	0.620	5.388	0.000	Frequency	ESQ-R total	Mean	0
-3.227	0.553	-5.831	0.000	Ability	ESQ-R total	Mean	0
6.701	1.067	6.282	0.000	Frequency	ESQ-R total	+1 SD	1
-4.923	0.873	-5.637	0.000	Ability	ESQ-R total	+1 SD	1
1.402	0.670	2.091	0.036	Frequency	EMQ-R total	-1 SD	-1
-1.372	0.666	-2.059	0.040	Ability	EMQ-R total	-1 SD	-1
2.956	0.521	5.678	0.000	Frequency	EMQ-R total	Mean	0
-1.694	0.426	-3.979	0.000	Ability	EMQ-R total	Mean	0
4.510	0.927	4.863	0.000	Frequency	EMQ-R total	+1 SD	1
-2.016	0.572	-3.527	0.000	Ability	EMQ-R total	+1 SD	1
0.428	1.157	0.370	0.711	Frequency	PRMQ total	-1 SD	-1
-0.212	1.500	-0.141	0.887	Ability	PRMQ total	-1 SD	-1
2.982	0.775	3.848	0.000	Frequency	PRMQ total	Mean	0
-1.807	0.907	-1.991	0.046	Ability	PRMQ total	Mean	0
5.535	0.941	5.883	0.000	Frequency	PRMQ total	+1 SD	1
-3.401	0.820	-4.150	0.000	Ability	PRMQ total	+1 SD	1

7 Interaction Figures

plot_data <- simple_slopes_agez %>%
  select(Outcome, Predictor, Age_probe, Age_z, b)

fig_freq <- plot_data %>%
  filter(Predictor == "Frequency") %>%
  ggplot(aes(x = Age_z, y = b, color = Outcome, group = Outcome)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_line(linewidth = 1) +
  geom_point(size = 2.5) +
  scale_x_continuous(breaks = c(-1, 0, 1), labels = c("-1 SD", "Mean", "+1 SD")) +
  labs(
    x = "Continuous age probe",
    y = "Frequency slope",
    title = "Age moderation of self-reported gaming frequency associations"
  ) +
  theme_minimal(base_size = 12)

fig_ability <- plot_data %>%
  filter(Predictor == "Ability") %>%
  ggplot(aes(x = Age_z, y = b, color = Outcome, group = Outcome)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_line(linewidth = 1) +
  geom_point(size = 2.5) +
  scale_x_continuous(breaks = c(-1, 0, 1), labels = c("-1 SD", "Mean", "+1 SD")) +
  labs(
    x = "Continuous age probe",
    y = "Perceived ability slope",
    title = "Age moderation of perceived gaming ability associations"
  ) +
  theme_minimal(base_size = 12)

ggsave("outputs/figures/Fig_Interactions_AgeFreq_R.png", fig_freq, width = 8, height = 5, dpi = 300)
ggsave("outputs/figures/Fig_Interactions_AgeAbility_R.png", fig_ability, width = 8, height = 5, dpi = 300)

fig_freq

fig_ability

8 Bonferroni-Corrected Correlations

primary_predictors <- c(
  "GSQ_Freq_Total_Rebuilt",
  "GSQ_Ability_Total_Rebuilt",
  "GSQ_Total_Skill_Rebuilt",
  freq_cols,
  ability_cols,
  legacy_skill_cols
)

secondary_outcomes <- c(
  "Plan_management",
  "Time_management",
  "Organization",
  "Emotional_regulation",
  "Behavioral_regulation"
)

correlation_rows <- function(data, predictors, outcomes, group_label) {
  expand_grid(predictor = predictors, outcome = names(outcomes)) %>%
    mutate(
      group = group_label,
      outcome_label = outcomes[outcome],
      n = map2_int(predictor, outcome, ~sum(complete.cases(data[[.x]], data[[.y]]))),
      r = map2_dbl(predictor, outcome, ~cor(data[[.x]], data[[.y]], use = "complete.obs")),
      p = map2_dbl(predictor, outcome, ~cor.test(data[[.x]], data[[.y]])$p.value)
    ) %>%
    select(group, predictor, outcome, outcome_label, n, r, p)
}

add_adjustments <- function(df) {
  df %>%
    mutate(
      p_bonferroni = p.adjust(p, method = "bonferroni"),
      p_holm = p.adjust(p, method = "holm"),
      p_fdr = p.adjust(p, method = "BH")
    )
}

overall_primary <- correlation_rows(dat, primary_predictors, primary_outcomes, "Overall") %>%
  add_adjustments()

agegroup_primary <- dat %>%
  group_split(AgeGroup) %>%
  map_dfr(function(d) {
    correlation_rows(d, primary_predictors, primary_outcomes, unique(d$AgeGroup))
  }) %>%
  add_adjustments()

all_outcomes <- c(primary_outcomes, setNames(secondary_outcomes, secondary_outcomes))

agegroup_all_outcomes <- dat %>%
  group_split(AgeGroup) %>%
  map_dfr(function(d) {
    correlation_rows(d, primary_predictors, all_outcomes, unique(d$AgeGroup))
  }) %>%
  add_adjustments()

legacy_predictors <- c(legacy_skill_cols, "GSQ_Total_Skill_Rebuilt")

legacy_age_x_freq <- dat %>%
  group_split(AgeGroup, Legacy_Frequency_Group) %>%
  map_dfr(function(d) {
    group_label <- paste0(
      unique(d$AgeGroup),
      "; legacy ",
      unique(d$Legacy_Frequency_Group),
      " frequency"
    )
    correlation_rows(d, legacy_predictors, all_outcomes, group_label)
  }) %>%
  add_adjustments()

wb <- openxlsx::createWorkbook()
openxlsx::addWorksheet(wb, "Overall_primary")
openxlsx::writeData(wb, "Overall_primary", overall_primary)
openxlsx::addWorksheet(wb, "AgeGroup_primary")
openxlsx::writeData(wb, "AgeGroup_primary", agegroup_primary)
openxlsx::addWorksheet(wb, "AgeGroup_all_outcomes")
openxlsx::writeData(wb, "AgeGroup_all_outcomes", agegroup_all_outcomes)
openxlsx::addWorksheet(wb, "Legacy_age_x_freq")
openxlsx::writeData(wb, "Legacy_age_x_freq", legacy_age_x_freq)
openxlsx::saveWorkbook(wb, "outputs/tables/Supp_Correlations_Bonferroni_R.xlsx", overwrite = TRUE)

overall_primary %>%
  filter(p_bonferroni < .05) %>%
  arrange(p_bonferroni) %>%
  kable(digits = 3)

group	predictor	outcome	outcome_label	n	r	p	p_bonferroni	p_holm	p_fdr
Overall	GSQ_Freq_Total_Rebuilt	EMQ-R	EMQ-R total	507	0.239	0.000	0.000	0.000	0.000
Overall	SPG_Freq	EMQ-R	EMQ-R total	507	0.207	0.000	0.000	0.000	0.000
Overall	GSQ_Total_Skill_Rebuilt	EMQ-R	EMQ-R total	507	0.187	0.000	0.001	0.001	0.000
Overall	STG_Freq	Total_Score_ESQ	ESQ-R total	507	0.184	0.000	0.002	0.002	0.000
Overall	SPG_CurrentSkill	EMQ-R	EMQ-R total	507	0.183	0.000	0.002	0.002	0.000
Overall	GSQ_Freq_Total_Rebuilt	PRMQ	PRMQ total	507	0.180	0.000	0.003	0.003	0.000
Overall	SPG_Freq	PRMQ	PRMQ total	507	0.179	0.000	0.003	0.003	0.000
Overall	AAG_Freq	Total_Score_ESQ	ESQ-R total	507	0.175	0.000	0.005	0.004	0.001
Overall	GSQ_Freq_Total_Rebuilt	Total_Score_ESQ	ESQ-R total	507	0.163	0.000	0.014	0.013	0.002
Overall	STG_Freq	EMQ-R	EMQ-R total	507	0.161	0.000	0.017	0.014	0.002
Overall	AAG_Freq	EMQ-R	EMQ-R total	507	0.161	0.000	0.018	0.015	0.002
Overall	FPSG_Freq	PRMQ	PRMQ total	507	0.152	0.001	0.037	0.030	0.003
Overall	SPG_CurrentSkill	PRMQ	PRMQ total	507	0.150	0.001	0.043	0.035	0.003

legacy_age_x_freq %>%
  filter(p_bonferroni < .05) %>%
  arrange(group, p_bonferroni) %>%
  kable(digits = 3)

group	predictor	outcome	outcome_label	n	r	p_bonferroni	p_holm	p_fdr
18-35; legacy High frequency	PSG_CurrentSkill	PRMQ	PRMQ total	103	0.419	0.003	0.003	0.001
18-35; legacy High frequency	RPG_CurrentSkill	Behavioral_regulation	Behavioral_regulation	103	-0.405	0.007	0.007	0.001
18-35; legacy High frequency	PSG_CurrentSkill	EMQ-R	EMQ-R total	103	0.387	0.018	0.017	0.002
36-55; legacy High frequency	AAG_CurrentSkill	Emotional_regulation	Emotional_regulation	123	-0.411	0.001	0.001	0.000
36-55; legacy High frequency	GSQ_Total_Skill_Rebuilt	Emotional_regulation	Emotional_regulation	123	-0.385	0.004	0.004	0.001
36-55; legacy High frequency	SPG_CurrentSkill	Emotional_regulation	Emotional_regulation	123	-0.350	0.025	0.024	0.003
36-55; legacy Low frequency	FPSG_CurrentSkill	Organization	Organization	95	-0.459	0.001	0.001	0.000
56+; legacy High frequency	SPG_CurrentSkill	PRMQ	PRMQ total	61	0.563	0.001	0.001	0.000
56+; legacy Low frequency	AAG_CurrentSkill	PRMQ	PRMQ total	42	-0.666	0.000	0.000	0.000
56+; legacy Low frequency	RPG_CurrentSkill	Organization	Organization	42	-0.562	0.036	0.035	0.004

9 Robustness Checks

Quadratic frequency terms are exploratory and should not replace the prespecified primary models.

The analyses below also reproduce the retained legacy/supplemental checks:

legacy summed GSQ current-skill models;
age-group presentation models;
legacy age-group x frequency correlation tables, reported above in the supplemental workbook.

The submitted manuscript’s data-driven “best” mixed-effects models are not reproduced here because they were removed from the revised manuscript rather than retained as supplemental evidence.

nonlinear_terms <- c(primary_terms, "GSQ_Freq_Total_z_sq", "Age_z_x_GSQ_Freq_Total_z_sq")
legacy_terms <- c("Sex_Male", "Education_z", "Age_z", "GSQ_Total_Skill_Rebuilt_z", "Age_z_x_GSQ_Total_Skill_z")
agegroup_terms <- c(
  "Sex_Male",
  "Education_z",
  "AgeGroup",
  "GSQ_Freq_Total_Rebuilt_z",
  "GSQ_Ability_Total_Rebuilt_z",
  "AgeGroup:GSQ_Freq_Total_Rebuilt_z",
  "AgeGroup:GSQ_Ability_Total_Rebuilt_z"
)

robustness_fits <- purrr::imap_dfr(primary_outcomes, function(label, outcome) {
  primary_formula <- as.formula(paste0("`", outcome, "` ~ ", paste(primary_terms, collapse = " + ")))
  nonlinear_formula <- as.formula(paste0("`", outcome, "` ~ ", paste(nonlinear_terms, collapse = " + ")))
  legacy_formula <- as.formula(paste0("`", outcome, "` ~ ", paste(legacy_terms, collapse = " + ")))
  agegroup_formula <- as.formula(paste0("`", outcome, "` ~ ", paste(agegroup_terms, collapse = " + ")))

  primary <- lm(primary_formula, data = dat)
  nonlinear <- lm(nonlinear_formula, data = dat)
  legacy <- lm(legacy_formula, data = dat)
  agegroup <- lm(agegroup_formula, data = dat)

  tibble(
    Outcome = label,
    Model = c(
      "Primary continuous-age",
      "Exploratory quadratic frequency",
      "Legacy summed skill",
      "Age-group presentation"
    ),
    n = c(nobs(primary), nobs(nonlinear), nobs(legacy), nobs(agegroup)),
    R2 = c(
      summary(primary)$r.squared,
      summary(nonlinear)$r.squared,
      summary(legacy)$r.squared,
      summary(agegroup)$r.squared
    ),
    AIC = c(AIC(primary), AIC(nonlinear), AIC(legacy), AIC(agegroup)),
    Delta_AIC_vs_primary = c(
      0,
      AIC(nonlinear) - AIC(primary),
      AIC(legacy) - AIC(primary),
      AIC(agegroup) - AIC(primary)
    )
  )
})

readr::write_csv(robustness_fits, "outputs/tables/Robustness_Model_Fit_R.csv")
kable(robustness_fits, digits = 3)

Outcome	Model	n	R2	AIC	Delta_AIC_vs_primary
ESQ-R total	Primary continuous-age	507	0.235	3787.603	0.000
ESQ-R total	Exploratory quadratic frequency	507	0.245	3785.004	-2.600
ESQ-R total	Legacy summed skill	507	0.132	3847.655	60.052
ESQ-R total	Age-group presentation	507	0.202	3815.226	27.623
EMQ-R total	Primary continuous-age	507	0.186	3579.204	0.000
EMQ-R total	Exploratory quadratic frequency	507	0.207	3569.766	-9.438
EMQ-R total	Legacy summed skill	507	0.123	3612.862	33.657
EMQ-R total	Age-group presentation	507	0.154	3604.913	25.709
PRMQ total	Primary continuous-age	507	0.131	3820.558	0.000
PRMQ total	Exploratory quadratic frequency	507	0.137	3821.031	0.474
PRMQ total	Legacy summed skill	507	0.067	3852.476	31.919
PRMQ total	Age-group presentation	507	0.126	3829.323	8.765

legacy_model_coefficients <- purrr::imap_dfr(primary_outcomes, function(label, outcome) {
  fit_lm_hc3(dat, outcome, legacy_terms) %>%
    filter(term != "(Intercept)") %>%
    mutate(
      Outcome = label,
      Model = "Legacy summed skill"
    )
})

agegroup_model_coefficients <- purrr::imap_dfr(primary_outcomes, function(label, outcome) {
  fit_lm_hc3(dat, outcome, agegroup_terms) %>%
    filter(term != "(Intercept)") %>%
    mutate(
      Outcome = label,
      Model = "Age-group presentation"
    )
})

retained_legacy_coefficients <- bind_rows(
  legacy_model_coefficients,
  agegroup_model_coefficients
) %>%
  mutate(
    p_holm_within_model_family = ave(p_raw, Model, FUN = function(p) p.adjust(p, method = "holm")),
    p_fdr_within_model_family = ave(p_raw, Model, FUN = function(p) p.adjust(p, method = "BH"))
  )

readr::write_csv(
  retained_legacy_coefficients,
  "outputs/tables/Retained_Legacy_And_AgeGroup_Coefficients_R.csv"
)

retained_legacy_coefficients %>%
  transmute(
    Outcome,
    Model,
    term,
    b = round(b, 3),
    SE_HC3 = round(robust_se, 3),
    p_raw = signif(p_raw, 3),
    p_holm_within_model_family = signif(p_holm_within_model_family, 3)
  ) %>%
  kable()

Outcome	Model	term	b	SE_HC3	p_raw	p_holm_within_model_family
ESQ-R total	Legacy summed skill	Sex_Male	0.258	1.005	7.98e-01	1.00e+00
ESQ-R total	Legacy summed skill	Education_z	-2.180	0.443	1.20e-06	1.77e-05
ESQ-R total	Legacy summed skill	Age_z	-2.786	0.681	5.04e-05	6.56e-04
ESQ-R total	Legacy summed skill	GSQ_Total_Skill_Rebuilt_z	0.244	0.620	6.95e-01	1.00e+00
ESQ-R total	Legacy summed skill	Age_z_x_GSQ_Total_Skill_z	1.810	0.953	5.81e-02	5.46e-01
EMQ-R total	Legacy summed skill	Sex_Male	0.284	0.858	7.41e-01	1.00e+00
EMQ-R total	Legacy summed skill	Education_z	-0.770	0.403	5.69e-02	5.46e-01
EMQ-R total	Legacy summed skill	Age_z	-2.161	0.515	3.24e-05	4.54e-04
EMQ-R total	Legacy summed skill	GSQ_Total_Skill_Rebuilt_z	1.204	0.533	2.44e-02	2.68e-01
EMQ-R total	Legacy summed skill	Age_z_x_GSQ_Total_Skill_z	1.180	0.722	1.03e-01	6.18e-01
PRMQ total	Legacy summed skill	Sex_Male	-1.176	1.109	2.89e-01	1.00e+00
PRMQ total	Legacy summed skill	Education_z	-0.939	0.487	5.46e-02	5.46e-01
PRMQ total	Legacy summed skill	Age_z	-1.913	0.598	1.46e-03	1.75e-02
PRMQ total	Legacy summed skill	GSQ_Total_Skill_Rebuilt_z	1.137	0.593	5.57e-02	5.46e-01
PRMQ total	Legacy summed skill	Age_z_x_GSQ_Total_Skill_z	1.133	0.800	1.57e-01	7.86e-01
ESQ-R total	Age-group presentation	Sex_Male	0.294	0.959	7.60e-01	1.00e+00
ESQ-R total	Age-group presentation	Education_z	-2.068	0.473	1.49e-05	3.73e-04
ESQ-R total	Age-group presentation	AgeGroup36-55	-5.542	1.025	1.00e-07	2.80e-06
ESQ-R total	Age-group presentation	AgeGroup56+	-11.300	1.490	0.00e+00	0.00e+00
ESQ-R total	Age-group presentation	GSQ_Freq_Total_Rebuilt_z	1.502	0.674	2.63e-02	4.73e-01
ESQ-R total	Age-group presentation	GSQ_Ability_Total_Rebuilt_z	-2.889	0.792	2.92e-04	6.71e-03
ESQ-R total	Age-group presentation	AgeGroup36-55:GSQ_Freq_Total_Rebuilt_z	0.185	1.343	8.90e-01	1.00e+00
ESQ-R total	Age-group presentation	AgeGroup56+:GSQ_Freq_Total_Rebuilt_z	6.230	1.935	1.37e-03	2.88e-02
ESQ-R total	Age-group presentation	AgeGroup36-55:GSQ_Ability_Total_Rebuilt_z	1.341	1.188	2.60e-01	1.00e+00
ESQ-R total	Age-group presentation	AgeGroup56+:GSQ_Ability_Total_Rebuilt_z	-2.629	1.526	8.56e-02	1.00e+00
EMQ-R total	Age-group presentation	Sex_Male	0.537	0.828	5.17e-01	1.00e+00
EMQ-R total	Age-group presentation	Education_z	-0.553	0.409	1.77e-01	1.00e+00
EMQ-R total	Age-group presentation	AgeGroup36-55	-3.940	0.878	8.90e-06	2.31e-04
EMQ-R total	Age-group presentation	AgeGroup56+	-6.782	1.243	1.00e-07	2.20e-06
EMQ-R total	Age-group presentation	GSQ_Freq_Total_Rebuilt_z	2.215	0.610	3.10e-04	6.81e-03
EMQ-R total	Age-group presentation	GSQ_Ability_Total_Rebuilt_z	-1.307	0.727	7.30e-02	1.00e+00
EMQ-R total	Age-group presentation	AgeGroup36-55:GSQ_Freq_Total_Rebuilt_z	0.054	1.076	9.60e-01	1.00e+00
EMQ-R total	Age-group presentation	AgeGroup56+:GSQ_Freq_Total_Rebuilt_z	1.818	1.588	2.53e-01	1.00e+00
EMQ-R total	Age-group presentation	AgeGroup36-55:GSQ_Ability_Total_Rebuilt_z	-0.149	0.969	8.78e-01	1.00e+00
EMQ-R total	Age-group presentation	AgeGroup56+:GSQ_Ability_Total_Rebuilt_z	0.281	1.083	7.96e-01	1.00e+00
PRMQ total	Age-group presentation	Sex_Male	-0.939	1.076	3.83e-01	1.00e+00
PRMQ total	Age-group presentation	Education_z	-0.812	0.468	8.33e-02	1.00e+00
PRMQ total	Age-group presentation	AgeGroup36-55	-4.243	1.060	7.21e-05	1.73e-03
PRMQ total	Age-group presentation	AgeGroup56+	-7.636	1.414	1.00e-07	2.80e-06
PRMQ total	Age-group presentation	GSQ_Freq_Total_Rebuilt_z	0.664	1.419	6.40e-01	1.00e+00
PRMQ total	Age-group presentation	GSQ_Ability_Total_Rebuilt_z	0.537	1.778	7.63e-01	1.00e+00
PRMQ total	Age-group presentation	AgeGroup36-55:GSQ_Freq_Total_Rebuilt_z	1.890	1.705	2.68e-01	1.00e+00
PRMQ total	Age-group presentation	AgeGroup56+:GSQ_Freq_Total_Rebuilt_z	6.285	2.008	1.85e-03	3.71e-02
PRMQ total	Age-group presentation	AgeGroup36-55:GSQ_Ability_Total_Rebuilt_z	-2.498	1.937	1.98e-01	1.00e+00
PRMQ total	Age-group presentation	AgeGroup56+:GSQ_Ability_Total_Rebuilt_z	-5.353	2.131	1.23e-02	2.34e-01

10 Session Information

sessionInfo()

## R version 4.4.1 (2024-06-14 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 11 x64 (build 26200)
## 
## Matrix products: default
## 
## 
## locale:
## [1] LC_COLLATE=English_United States.utf8 
## [2] LC_CTYPE=English_United States.utf8   
## [3] LC_MONETARY=English_United States.utf8
## [4] LC_NUMERIC=C                          
## [5] LC_TIME=English_United States.utf8    
## 
## time zone: Europe/Athens
## tzcode source: internal
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] knitr_1.50      openxlsx_4.2.8  lmtest_0.9-40   zoo_1.8-12     
##  [5] sandwich_3.1-1  broom_1.0.8     readxl_1.4.5    lubridate_1.9.4
##  [9] forcats_1.0.0   stringr_1.5.1   dplyr_1.1.4     purrr_1.0.4    
## [13] readr_2.1.5     tidyr_1.3.1     tibble_3.2.1    ggplot2_4.0.0  
## [17] tidyverse_2.0.0
## 
## loaded via a namespace (and not attached):
##  [1] sass_0.4.10        generics_0.1.4     stringi_1.8.4      lattice_0.22-6    
##  [5] hms_1.1.3          digest_0.6.37      magrittr_2.0.3     timechange_0.3.0  
##  [9] evaluate_1.0.4     grid_4.4.1         RColorBrewer_1.1-3 fastmap_1.2.0     
## [13] cellranger_1.1.0   jsonlite_2.0.0     zip_2.3.3          backports_1.5.0   
## [17] scales_1.4.0       textshaping_1.0.1  jquerylib_0.1.4    cli_3.6.3         
## [21] crayon_1.5.3       rlang_1.1.4        bit64_4.6.0-1      withr_3.0.2       
## [25] cachem_1.1.0       yaml_2.3.10        parallel_4.4.1     tools_4.4.1       
## [29] tzdb_0.5.0         vctrs_0.6.5        R6_2.6.1           lifecycle_1.0.4   
## [33] bit_4.6.0          vroom_1.6.5        ragg_1.4.0         pkgconfig_2.0.3   
## [37] pillar_1.10.2      bslib_0.9.0        gtable_0.3.6       glue_1.7.0        
## [41] Rcpp_1.0.13        systemfonts_1.2.3  xfun_0.52          tidyselect_1.2.1  
## [45] rstudioapi_0.17.1  farver_2.1.2       htmltools_0.5.8.1  labeling_0.4.3    
## [49] rmarkdown_2.29     compiler_4.4.1     S7_0.2.0

Analyses for Moderation is Best Paper

Gaming frequency, perceived gaming ability, and self-reported everyday cognition across adulthood

Lepida, Nega, MacPherson, & Kourtesis

2026-05-04