# Maximal respiration (FCCP minus rotenone/antimycin A) in KO vs WT: raw OCR with wells as n, OCR per 1,000 cells with
# wells as n, and OCR per 1,000 cells with plates as n. R 4.3.3, base R.
d <- read.csv("seahorse.csv")                        # one row per well: plate, genotype, well, cells, ocr_basal, ocr_oligomycin, ocr_fccp, ocr_rot_aa
d$maximal <- d$ocr_fccp - d$ocr_rot_aa; d$maximal_per_k <- d$maximal / d$cells * 1000
ko <- d$genotype == "KO"; f2 <- function(p) formatC(signif(p, 2), format = "g")
wells <- function(col, lab) {
  t <- t.test(log(d[[col]][ko]), log(d[[col]][!ko]))
  cat(sprintf("%s, wells as n (%d vs %d): KO / WT %.2f, p = %s\n", lab, sum(ko), sum(!ko), exp(mean(log(d[[col]][ko])) - mean(log(d[[col]][!ko]))), f2(t$p.value)))
}
wells("maximal", "Maximal respiration, raw OCR")
wells("maximal_per_k", "Maximal respiration per 1,000 cells")
plates <- function(col, lab) {
  m <- tapply(log(d[[col]]), list(d$plate, d$genotype), mean); r <- m[, "KO"] - m[, "WT"]; t <- t.test(r)
  cat(sprintf("%s, plates as n (%d plates): KO / WT %.2f (95%% CI %.2f to %.2f), p = %s; per-plate ratios %s\n", lab, length(r), exp(mean(r)),
              exp(t$conf.int[1]), exp(t$conf.int[2]), f2(t$p.value), paste(sprintf("%.2f", exp(r)), collapse = ", ")))
}
plates("maximal", "Maximal respiration, raw OCR")
plates("maximal_per_k", "Maximal respiration per 1,000 cells")
cat(sprintf("Mean cells per well: WT %.0f, KO %.0f\n", mean(d$cells[!ko]), mean(d$cells[ko])))
