# Aggregate results over the pre-registered associations. Proportions carry exact
# (Clopper-Pearson) 95% intervals and medians the distribution-free interval from order
# statistics. Values are rounded to the precision their uncertainty supports, so a re-run on
# another platform declares the same numbers.

share <- function(k, n) {
  if (n == 0) return(list(k = k, n = n, share = NA, ci = c(NA, NA)))
  list(k = k, n = n, share = round(k / n, 3), ci = round(as.numeric(stats::binom.test(k, n)$conf.int), 3))
}

median_with_ci <- function(x) {
  x <- x[is.finite(x)]
  if (length(x) == 0) return(list(n = 0, median = NA, ci = c(NA, NA), quartiles = c(NA, NA)))
  list(n = length(x), median = round(stats::median(x), 3), ci = round(median_ci(x), 3),
       quartiles = round(as.numeric(stats::quantile(x, c(0.25, 0.75), type = 7)), 3))
}

# The departures the coding check found bearing on the headline estimate: how many associations
# have one, by kind, and how many have one the re-implementation doesn't follow (it corrects it,
# or the paper's numbers rule out the stated computation without revealing its own).
departure_counts <- function(results) {
  n <- length(results)
  count <- function(f) sum(vapply(results, function(r) isTRUE(f(r)), TRUE))
  affecting <- function(r, keep = function(d) TRUE) {
    any(vapply(r$departures, function(d) isTRUE(d$affects_headline) && keep(d), TRUE))
  }
  list(
    affecting_headline = share(count(function(r) affecting(r)), n),
    by_kind = setNames(lapply(DEPARTURE_KINDS, function(k) share(count(function(r) affecting(r, function(d) d$kind == k)), n)), DEPARTURE_KINDS),
    not_followed = share(count(function(r) affecting(r, function(d) isFALSE(d$followed))), n)
  )
}

summarize_results <- function(results, replication) {
  count <- function(rs, f) sum(vapply(rs, function(r) isTRUE(f(r)), TRUE))
  n <- length(results)
  out <- list(reproduction = list(
    associations = n,
    reproduced = share(count(results, function(r) r$reproduced), n),
    same_sign_p05 = share(count(results, function(r) r$same_sign_original && r$original$p < 0.05), n),
    ratio_original = median_with_ci(vapply(results, function(r) r$ratio_original, 0)),
    departures = departure_counts(results)
  ))
  if (replication) {
    informative <- Filter(function(r) isTRUE(r$informative), results)
    reproduced <- Filter(function(r) isTRUE(r$reproduced), results)
    out$replication <- list(
      associations = n,
      replicated = share(count(results, function(r) r$replicated), n),
      informative = length(informative),
      replicated_informative = share(count(informative, function(r) r$replicated), length(informative)),
      same_sign_p05 = share(count(results, function(r) r$same_sign && r$replication$p < 0.05), n),
      in_published_ci = share(count(results, function(r) r$in_published_ci), n),
      reversed = share(count(results, function(r) !r$same_sign && r$q < 0.05), n),
      ratio = median_with_ci(vapply(results, function(r) r$ratio, 0)),
      ratio_informative = median_with_ci(vapply(informative, function(r) r$ratio, 0)),
      ratio_own = median_with_ci(vapply(results, function(r) r$ratio_own, 0)),
      heterogeneous = share(count(results, function(r) r$heterogeneity_q < 0.05), n),
      differs_from_published = share(count(results, function(r) r$differs_from_published), n),
      differs_by_direction = lapply(c(smaller = "smaller", larger = "larger", opposite_sign = "opposite sign"),
                                    function(d) share(count(results, function(r) identical(r$difference, d)), n)),
      replicated_reproduced = share(count(reproduced, function(r) r$replicated), length(reproduced)),
      replicated_paper_weight = share(count(results, function(r) r$replicated_paper_weight), n)
    )
  }
  out
}
