Correlation Analysis — What Moves Together

suppressPackageStartupMessages(library(DT))
suppressPackageStartupMessages(library(htmlwidgets))
suppressPackageStartupMessages(library(arrow))
suppressPackageStartupMessages(library(knitr))
suppressPackageStartupMessages(library(rmarkdown))
suppressPackageStartupMessages(library(dplyr))
suppressPackageStartupMessages(library(tidyr))
suppressPackageStartupMessages(library(ggplot2))
suppressPackageStartupMessages(library(stringr))
suppressPackageStartupMessages(library(lubridate))
suppressPackageStartupMessages(library(broom))
suppressPackageStartupMessages(library(Matrix))
suppressPackageStartupMessages(library(cluster))
suppressPackageStartupMessages(library(data.table))

compute_shared <- function(df, params, col_map = list()) {
  # === SHARED EXPORTS ===
  #   initial_rows/final_rows/rows_removed  $ row accounting
  #   feature_names     $ named character — semantic -> humanized names
  #   used_features     $ character — semantic names used
  #   dropped_features  $ character — excluded columns
  #   cor_long_df       $ data.frame(feature_a, feature_b, correlation) — full matrix
  #   pairs_df          $ data.frame(pair, correlation, p_value, n, strength, significance)
  #   top_pair_df       $ data.frame(value_a, value_b) — <=1000 sample of strongest pair
  #   top_pair_names    $ character(2) — humanized names of the strongest pair
  #   top_r / top_p     $ numeric — strongest pair stats
  #   n_sig             $ integer — pairs significant at p<0.05
  #   metrics / json_output
  # === /SHARED EXPORTS ===

initial_rows <- nrow(df)
  feat_cols <- grep("^feature_[0-9]+$", names(df), value = TRUE)
  feat_cols <- feat_cols[order(as.integer(sub("^feature_", "", feat_cols)))]
  if (length(feat_cols) < 2) {
    stop("column_mapping must map at least two feature columns(feature_1, feature_2)")
  }
  feature_names <- setNames(humanize_semantic(feat_cols, col_map), feat_cols)

dropped_features <- character(0)
  for (fc in feat_cols) {
    v <- df[[fc]]
    if (!is.numeric(v)) {
      conv <- suppressWarnings(as.numeric(as.character(v)))
      n_orig <- sum(!is.na(v) & as.character(v) != "")
      if (n_orig > 0 && sum(!is.na(conv)) >= 0.95 * n_orig) {
        df[[fc]] <- conv
      } else {
        dropped_features <- c(dropped_features, fc); next
      }
    }
    v <- df[[fc]]
    med <- median(v, na.rm = TRUE)
    if (is.na(med)) { dropped_features <- c(dropped_features, fc); next }
    v[is.na(v)] <- med
    df[[fc]] <- v
    if (isTRUE(var(v) == 0) || is.na(var(v))) {
      dropped_features <- c(dropped_features, fc)
    }
  }
  used_features <- setdiff(feat_cols, dropped_features)
  if (length(used_features) < 2) {
    stop(paste0("Correlation needs at least two usable numeric columns; only ",
                length(used_features), " remained after cleaning."))
  }
  X <- df[, used_features, drop = FALSE]
  final_rows <- nrow(X)
  rows_removed <- initial_rows - final_rows
  if (final_rows < 10) stop(sprintf("Only %d usable rows — need at least 10.", final_rows))

C <- suppressWarnings(cor(X, use = "pairwise.complete.obs", method = "pearson"))
  hn <- unname(feature_names[used_features])

k <- length(used_features)
  cor_long_df <- do.call(rbind, lapply(seq_len(k), function(i) {
    data.frame(feature_a = hn[i], feature_b = hn,
               correlation = round(as.numeric(C[i, ]), 3),
               stringsAsFactors = FALSE)
  }))
  rownames(cor_long_df) <- NULL

pair_rows <- list()
  for (i in seq_len(k - 1)) {
    for (j in (i + 1):k) {
      r_ij <- C[i, j]
      if (is.na(r_ij)) next
      ct <- tryCatch(cor.test(X[[i]], X[[j]]), error = function(e) NULL)
      p_ij <- if (!is.null(ct)) ct$p.value else NA_real_
      pair_rows[[length(pair_rows) + 1]] <- data.frame(
        pair = paste0(hn[i], " ~ ", hn[j]),
        a = used_features[i], b = used_features[j],
        correlation = round(r_ij, 3),
        p_value = signif(p_ij, 3),
        n = sum(complete.cases(X[[i]], X[[j]])),
        stringsAsFactors = FALSE
      )
    }
  }
  pairs_all <- do.call(rbind, pair_rows)
  if (is.null(pairs_all) || nrow(pairs_all) == 0) {
    stop("No valid correlation pairs could be computed from the mapped columns.")
  }
  pairs_all <- pairs_all[order(-abs(pairs_all$correlation)), , drop = FALSE]
  rownames(pairs_all) <- NULL
  pairs_all$strength <- ifelse(abs(pairs_all$correlation) >= 0.7, "strong",
                        ifelse(abs(pairs_all$correlation) >= 0.4, "moderate", "weak"))
  pairs_all$significance <- ifelse(is.na(pairs_all$p_value), "",
                            ifelse(pairs_all$p_value < 0.001, "***",
                            ifelse(pairs_all$p_value < 0.01, "**",
                            ifelse(pairs_all$p_value < 0.05, "*", ""))))
  n_sig <- sum(pairs_all$p_value < 0.05, na.rm = TRUE)
  pairs_df <- head(pairs_all[, c("pair", "correlation", "p_value", "n",
                                 "strength", "significance")], 15)

top <- pairs_all[1, ]
  top_pair_names <- c(feature_names[[top$a]], feature_names[[top$b]])
  set.seed(42)
  sidx <- if (final_rows > 1000) sample(final_rows, 1000) else seq_len(final_rows)
  top_pair_df <- data.frame(
    value_a = X[[top$a]][sidx],
    value_b = X[[top$b]][sidx],
    stringsAsFactors = FALSE
  )
  top_pair_df <- top_pair_df[order(top_pair_df$value_a), ]
  rownames(top_pair_df) <- NULL
  top_r <- top$correlation
  top_p <- top$p_value

  metrics <- list(
    `Observations`        = final_rows,
    `Columns Analysed`    = length(used_features),
    `Pairs Tested`        = nrow(pairs_all),
    `Significant Pairs`   = as.integer(n_sig),
    `Strongest |r|`       = round(abs(top_r), 3),
    `Strongest Pair`      = top$pair
  )

  json_output <- list(
    answer = paste0(
      "Pairwise Pearson correlations on ", length(used_features),
      " columns across ", format(final_rows, big.mark = ","), " rows: ",
      "strongest pair is ", top$pair, " (r = ", top_r,
      if (!is.na(top_p)) paste0(", p = ", top_p) else "", "); ",
      n_sig, " of ", nrow(pairs_all), " pairs significant at p<0.05."
    ),
    cards = lapply(
      c("tldr", "overview", "preprocessing", "correlation_heatmap",
        "strongest_pairs", "top_pair_scatter"),
      function(cid) list(id = cid, metrics = metrics)
    )
  )

  list(
    initial_rows = initial_rows, final_rows = final_rows,
    rows_removed = rows_removed,
    feature_names = feature_names, used_features = used_features,
    dropped_features = dropped_features,
    cor_long_df = cor_long_df, pairs_df = pairs_df, pairs_all = pairs_all,
    top_pair_df = top_pair_df, top_pair_names = top_pair_names,
    top_r = top_r, top_p = top_p, n_sig = n_sig,
    metrics = metrics, json_output = json_output
  )
}