# 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. Python 3.13, NumPy, pandas, SciPy.
import numpy as np, pandas as pd
from scipy import stats

d = pd.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["basal"] = d.ocr_basal - d.ocr_rot_aa
d["maximal_per_k"] = d.maximal / d.cells * 1000; d["basal_per_k"] = d.basal / d.cells * 1000
ko = d.genotype == "KO"
def wells(col, lab):
    t = stats.ttest_ind(np.log(d[col][ko]), np.log(d[col][~ko]), equal_var=False)
    print(f"{lab}, wells as n ({ko.sum()} vs {(~ko).sum()}): KO / WT {np.exp(np.log(d[col][ko]).mean() - np.log(d[col][~ko]).mean()):.2f}, p = {t.pvalue:.2g}")
wells("maximal", "Maximal respiration, raw OCR")
wells("maximal_per_k", "Maximal respiration per 1,000 cells")
def plates(col, lab):
    m = d.groupby(["plate", "genotype"])[col].apply(lambda v: np.log(v).mean()).unstack()
    r = (m["KO"] - m["WT"]).to_numpy(); t = stats.ttest_1samp(r, 0); ci = t.confidence_interval()
    print(f"{lab}, plates as n ({len(r)} plates): KO / WT {np.exp(r.mean()):.2f} (95% CI {np.exp(ci.low):.2f} to {np.exp(ci.high):.2f}), p = {t.pvalue:.2g}; "
          f"per-plate ratios " + ", ".join(f"{v:.2f}" for v in np.exp(r)))
plates("maximal", "Maximal respiration, raw OCR")
plates("maximal_per_k", "Maximal respiration per 1,000 cells")
print(f"Mean cells per well: WT {d.cells[~ko].mean():.0f}, KO {d.cells[ko].mean():.0f}")