---
title: "Effect size computation (Pinillos 2024 — Bank Cases, Stakes, and Normative Facts)"
format:
  html:
    toc: true
execute:
  echo: true
  warning: true
  message: false
---

## How to use

1) Copy this file once to `papers/<paper_key>/analysis/effect_sizes.qmd`.
2) For each computed effect, fill an **Inputs** block (duplicate the section and change `study_id`/`effect_id` as needed).
3) Render with Quarto (or run the R chunks in your IDE) and copy the resulting `d` and `v` into the paper YAML.
4) Keep the rendered HTML as a human-auditable record: it reports the method used + key inputs.

## Inputs

```{r}
# Identify the record you are computing.
paper_key <- "pinillos2024bankcasesstakes"
study_id <- 1
effect_id <- "s1_e1"

# Choose one method (make it explicit; this will be printed in the audit output).
#
# Between-subjects (high vs low are different participants):
# - "between_groups": means + SDs + n_high/n_low
# - "between_t": t + n_high/n_low (sign from means or sign_d)
# - "between_f": F + n_high/n_low (sign from means or sign_d)
# - "between_chi2": chi2 + n_high/n_low (sign from means or sign_d)
# - "between_2x2_or": 2x2 counts (yes/no by low/high) via esc::esc_2x2 (OR -> d)
# - "between_reported_d_t_df": reported d + t(df) (computes v without split Ns)
#
# Within-subjects (same participants rate both cases):
# - "within_smcrp_r": means + SDs + n_total + within-person r
# - "within_smcrp_t": means + SDs + paired t(df) to recover r
method_used <- "between_2x2_or"

# Sign convention (matches extraction instructions).
sign_convention <- "d = mean(low) - mean(high)"

# Required sample sizes (per condition).
n_high <- NA_integer_
n_low  <- NA_integer_

# For within-subject designs, use n_total instead (or provide df and we infer n_total = df+1).
n_total <- NA_integer_

# If available (recommended for sign): condition-level means and SDs
mean_high <- NA_real_
mean_low  <- NA_real_
sd_high   <- NA_real_
sd_low    <- NA_real_

# For binary outcomes, provide exact 2x2 counts when available.
yes_low  <- 35
no_low   <- 14
yes_high <- 5
no_high  <- 35

# If available: test statistics for 2-group contrasts (independent or paired)
t_value    <- NA_real_
f_value    <- NA_real_
chi2_value <- NA_real_
df         <- NA_real_

# If the paper reports d (sometimes without enough to compute v via group Ns), enter it here.
reported_d <- NA_real_

# For within-subject designs, provide r directly when known.
r_within <- NA_real_

# If the direction cannot be inferred from means, set this manually to +1 or -1.
# Convention: d = mean(low) - mean(high)
sign_d <- NA_real_

# Human-auditable metadata (printed in the audit output).
notes_on_assumptions <- "Binary DV (endorse 'Bob knows' vs not). Use exact 2x2 counts and esc::esc_2x2 (OR -> d), with group1=low stakes and group2=high stakes to match d = mean(low) - mean(high)."
imputed_flag <- FALSE
needs_sensitivity <- TRUE
```

## Computation

```{r}
stop_if_missing <- function(x, name) {
  if (is.na(x)) stop(sprintf("Missing required input: %s", name), call. = FALSE)
}

# Exact small-sample correction factor used by metafor (.cmicalc).
hedges_correction <- function(df) {
  ifelse(df <= 1, NA_real_, exp(lgamma(df/2) - log(sqrt(df/2)) - lgamma((df - 1)/2)))
}

pooled_sd <- function(n_high, n_low, sd_high, sd_low) {
  sqrt(((n_high - 1) * sd_high^2 + (n_low - 1) * sd_low^2) / (n_high + n_low - 2))
}

d_from_groups_independent <- function(n_high, n_low, mean_high, mean_low, sd_high, sd_low) {
  s <- pooled_sd(n_high, n_low, sd_high, sd_low)
  (mean_low - mean_high) / s
}

d_from_t_independent <- function(t_value, n_high, n_low) {
  t_value * sqrt((n_high + n_low) / (n_high * n_low))
}

d_from_f_independent <- function(f_value, n_high, n_low) {
  # Assumes a 2-group contrast where F = t^2
  sqrt(f_value * (n_high + n_low) / (n_high * n_low))
}

d_from_chi2_independent <- function(chi2_value, n_high, n_low) {
  sqrt(chi2_value * (n_high + n_low) / (n_high * n_low))
}

var_d_independent <- function(d, n_high, n_low) {
  n <- n_high + n_low
  (n / (n_high * n_low)) + (d^2 / (2 * (n - 2)))
}

# Between-subjects: compute v without split group Ns using reported d + t(df).
var_d_between_from_d_t_df <- function(d, t_value, df) {
  (d / t_value)^2 + (d^2) / (2 * df)
}

# Within-subjects (SMCRP-style): standardize by pooled SD across occasions.
sd_pooled_within <- function(sd_low, sd_high) {
  sqrt((sd_low^2 + sd_high^2) / 2)
}

d_within_smcrp <- function(mean_low, mean_high, sd_low, sd_high) {
  (mean_low - mean_high) / sd_pooled_within(sd_low, sd_high)
}

# Recover within-person correlation r from paired t + means/SDs (when r is not reported).
r_from_paired_t <- function(mean_low, mean_high, sd_low, sd_high, t_value, n_total) {
  mean_diff <- mean_low - mean_high
  sd_diff <- abs(mean_diff) * sqrt(n_total) / abs(t_value)
  (sd_low^2 + sd_high^2 - sd_diff^2) / (2 * sd_low * sd_high)
}

var_d_within_smcrp <- function(d, r, n_total) {
  (2 * (1 - r) / n_total) + (d^2) * (1 + r^2) / (4 * n_total)
}

infer_sign <- function(mean_low, mean_high, sign_d) {
  if (!is.na(mean_low) && !is.na(mean_high) && mean_low != mean_high) {
    return(sign(mean_low - mean_high))
  }
  stop_if_missing(sign_d, "sign_d (cannot infer sign from means)")
  if (!is.element(sign_d, c(-1, 1))) stop("sign_d must be +1 or -1", call. = FALSE)
  sign_d
}

# Compute d/v (and small-sample corrected g/v_g) based on method_used.
d <- NA_real_
v <- NA_real_
g <- NA_real_
v_g <- NA_real_
computed_from_suggested <- NA_character_
design_used <- if (startsWith(method_used, "between_")) "Between-Subjects" else if (startsWith(method_used, "within_")) "Within-Subjects" else NA_character_

if (method_used == "between_groups") {
  computed_from_suggested <- "groups"
  stop_if_missing(n_high, "n_high")
  stop_if_missing(n_low, "n_low")
  stop_if_missing(mean_high, "mean_high")
  stop_if_missing(mean_low, "mean_low")
  stop_if_missing(sd_high, "sd_high")
  stop_if_missing(sd_low, "sd_low")
  d <- d_from_groups_independent(n_high, n_low, mean_high, mean_low, sd_high, sd_low)
  v <- var_d_independent(d, n_high, n_low)
  df_used <- n_high + n_low - 2
  J <- hedges_correction(df_used)
  g <- J * d
  v_g <- (J^2) * v
} else if (method_used == "between_2x2_or") {
  computed_from_suggested <- "groups"
  stop_if_missing(yes_low, "yes_low")
  stop_if_missing(no_low, "no_low")
  stop_if_missing(yes_high, "yes_high")
  stop_if_missing(no_high, "no_high")
  if (!requireNamespace("esc", quietly = TRUE)) {
    stop("Package 'esc' is required for method_used='between_2x2_or'.", call. = FALSE)
  }

  esc_fit <- esc::esc_2x2(
    grp1yes = yes_low,
    grp1no = no_low,
    grp2yes = yes_high,
    grp2no = no_high,
    es.type = "d"
  )

  d <- as.numeric(esc_fit$es)
  v <- as.numeric(esc_fit$var)
  n_low <- yes_low + no_low
  n_high <- yes_high + no_high
  mean_low <- yes_low / n_low
  mean_high <- yes_high / n_high
  df_used <- n_high + n_low - 2
  J <- hedges_correction(df_used)
  g <- J * d
  v_g <- (J^2) * v
} else if (method_used == "between_t") {
  computed_from_suggested <- "t_df"
  stop_if_missing(n_high, "n_high")
  stop_if_missing(n_low, "n_low")
  stop_if_missing(t_value, "t_value")
  sign_used <- infer_sign(mean_low, mean_high, sign_d)
  d <- sign_used * abs(d_from_t_independent(t_value, n_high, n_low))
  v <- var_d_independent(d, n_high, n_low)
  df_used <- n_high + n_low - 2
  J <- hedges_correction(df_used)
  g <- J * d
  v_g <- (J^2) * v
} else if (method_used == "between_f") {
  computed_from_suggested <- "f_df"
  stop_if_missing(n_high, "n_high")
  stop_if_missing(n_low, "n_low")
  stop_if_missing(f_value, "f_value")
  sign_used <- infer_sign(mean_low, mean_high, sign_d)
  d <- sign_used * abs(d_from_f_independent(f_value, n_high, n_low))
  v <- var_d_independent(d, n_high, n_low)
  df_used <- n_high + n_low - 2
  J <- hedges_correction(df_used)
  g <- J * d
  v_g <- (J^2) * v
} else if (method_used == "between_chi2") {
  computed_from_suggested <- "other"
  stop_if_missing(n_high, "n_high")
  stop_if_missing(n_low, "n_low")
  stop_if_missing(chi2_value, "chi2_value")
  sign_used <- infer_sign(mean_low, mean_high, sign_d)
  d <- sign_used * abs(d_from_chi2_independent(chi2_value, n_high, n_low))
  v <- var_d_independent(d, n_high, n_low)
  df_used <- n_high + n_low - 2
  J <- hedges_correction(df_used)
  g <- J * d
  v_g <- (J^2) * v
} else if (method_used == "between_reported_d_t_df") {
  computed_from_suggested <- "reported_d"
  stop_if_missing(reported_d, "reported_d")
  stop_if_missing(t_value, "t_value")
  stop_if_missing(df, "df")
  sign_used <- infer_sign(mean_low, mean_high, sign_d)
  d <- sign_used * abs(reported_d)
  v <- var_d_between_from_d_t_df(d = d, t_value = abs(t_value), df = df)
  J <- hedges_correction(df)
  g <- J * d
  v_g <- (J^2) * v
} else if (method_used == "within_smcrp_r") {
  computed_from_suggested <- "groups"
  stop_if_missing(n_total, "n_total")
  stop_if_missing(mean_high, "mean_high")
  stop_if_missing(mean_low, "mean_low")
  stop_if_missing(sd_high, "sd_high")
  stop_if_missing(sd_low, "sd_low")
  stop_if_missing(r_within, "r_within")
  if (abs(r_within) > 1) stop("r_within must be between -1 and 1", call. = FALSE)
  d <- d_within_smcrp(mean_low, mean_high, sd_low, sd_high)
  v <- var_d_within_smcrp(d = d, r = r_within, n_total = n_total)
  df_used <- 2 * (n_total - 1) / (1 + r_within^2)
  J <- hedges_correction(df_used)
  g <- J * d
  v_g <- (2 * (1 - r_within) / n_total) + (g^2) * (1 + r_within^2) / (4 * n_total)
} else if (method_used == "within_smcrp_t") {
  computed_from_suggested <- "groups"
  stop_if_missing(mean_high, "mean_high")
  stop_if_missing(mean_low, "mean_low")
  stop_if_missing(sd_high, "sd_high")
  stop_if_missing(sd_low, "sd_low")
  stop_if_missing(t_value, "t_value")

  if (is.na(n_total)) {
    stop_if_missing(df, "df (or provide n_total)")
    n_total <- df + 1
  }
  if (!is.na(df) && df != (n_total - 1)) {
    warning("For a paired t-test, df should equal n_total - 1; check inputs.")
  }

  r_est <- r_from_paired_t(mean_low, mean_high, sd_low, sd_high, t_value, n_total)
  if (abs(r_est) > 1) stop(sprintf("Recovered r=%.4f outside [-1,1]; check inputs or provide r_within directly.", r_est), call. = FALSE)

  d <- d_within_smcrp(mean_low, mean_high, sd_low, sd_high)
  v <- var_d_within_smcrp(d = d, r = r_est, n_total = n_total)
  df_used <- 2 * (n_total - 1) / (1 + r_est^2)
  J <- hedges_correction(df_used)
  g <- J * d
  v_g <- (2 * (1 - r_est) / n_total) + (g^2) * (1 + r_est^2) / (4 * n_total)
} else {
  stop(sprintf("Unknown method_used: %s", method_used), call. = FALSE)
}

# Build a compact, human-auditable summary of inputs actually used.
inputs_used <- paste(
  c(
    sprintf("method=%s", method_used),
    sprintf("sign_convention=%s", sign_convention),
    if (!is.na(n_low)) sprintf("n_low=%s", n_low) else NULL,
    if (!is.na(n_high)) sprintf("n_high=%s", n_high) else NULL,
    if (!is.na(n_total)) sprintf("n_total=%s", n_total) else NULL,
    if (!is.na(df)) sprintf("df=%s", df) else NULL,
    if (!is.na(mean_low)) sprintf("mean_low=%s", mean_low) else NULL,
    if (!is.na(mean_high)) sprintf("mean_high=%s", mean_high) else NULL,
    if (!is.na(sd_low)) sprintf("sd_low=%s", sd_low) else NULL,
    if (!is.na(sd_high)) sprintf("sd_high=%s", sd_high) else NULL,
    if (!is.na(yes_low)) sprintf("yes_low=%s", yes_low) else NULL,
    if (!is.na(no_low)) sprintf("no_low=%s", no_low) else NULL,
    if (!is.na(yes_high)) sprintf("yes_high=%s", yes_high) else NULL,
    if (!is.na(no_high)) sprintf("no_high=%s", no_high) else NULL,
    if (!is.na(t_value)) sprintf("t=%s", t_value) else NULL,
    if (!is.na(f_value)) sprintf("f=%s", f_value) else NULL,
    if (!is.na(chi2_value)) sprintf("chi2=%s", chi2_value) else NULL,
    if (!is.na(reported_d)) sprintf("reported_d=%s", reported_d) else NULL,
    if (!is.na(r_within)) sprintf("r_within=%s", r_within) else NULL
  ),
  collapse = ", "
)

# Audit row (recommended to keep in rendered HTML).
audit <- data.frame(
  paper_key = paper_key,
  study_id = study_id,
  effect_id = effect_id,
  design = design_used,
  method_used = method_used,
  computed_from_suggested = computed_from_suggested,
  inputs_used = inputs_used,
  d = d,
  v = v,
  g = g,
  v_g = v_g,
  notes_on_assumptions = notes_on_assumptions,
  imputed_flag = imputed_flag,
  needs_sensitivity = needs_sensitivity
)
audit
```

## Paste-ready YAML snippet

```{r}
cat(sprintf(
  "effect_size:\\n  metric: SMD\\n  d: %.12f\\n  v: %.12f\\n  computed_from: %s\\n  needs_review: false\\n  notes: \"%s\"\\n",
  d, v, computed_from_suggested, gsub("\"", "'", inputs_used)
))
```