# Bliss independence and highest single agent (HSA) for a single-dose combination, per plate and across experiments.
# Python 3.13, NumPy, pandas, SciPy.
import numpy as np, pandas as pd
from scipy import stats

d = pd.read_csv("viability.csv")                     # columns: experiment, treatment (control, A, B, A+B), well, signal
d["sf"] = d.signal / d.groupby("experiment").signal.transform(lambda s: s[d.loc[s.index, "treatment"] == "control"].mean())

def plate(e):                                         # surviving fractions on one plate: mean well signal / mean control signal
    x = d[d.experiment == e]; return {g: x.sf[x.treatment == g].mean() for g in ("A", "B", "A+B")}, x

sf, x = plate(1)
bliss, hsa = sf["A"] * sf["B"], min(sf["A"], sf["B"])
print(f"Experiment 1: surviving fraction A {sf['A']:.3f}, B {sf['B']:.3f}, A+B {sf['A+B']:.3f}; Bliss expectation {bliss:.3f}; HSA (best single drug) {hsa:.3f}")
best = "A" if sf["A"] < sf["B"] else "B"
t1 = stats.ttest_ind(x.sf[x.treatment == "A+B"], x.sf[x.treatment == best], equal_var=False)
t2 = stats.ttest_1samp(x.sf[x.treatment == "A+B"], bliss)
print(f"  Wells as n, HSA: A+B vs {best} alone, Welch's t-test p = {t1.pvalue:.2g}")
print(f"  Wells as n, Bliss: A+B wells vs the Bliss expectation, one-sample t-test p = {t2.pvalue:.2g}")

print("Each experiment: observed / Bliss-expected surviving fraction [observed / HSA]:")
lb, lh = [], []
for e in sorted(d.experiment.unique()):
    sf, _ = plate(e)
    lb.append(np.log(sf["A+B"]) - np.log(sf["A"]) - np.log(sf["B"])); lh.append(np.log(sf["A+B"]) - np.log(min(sf["A"], sf["B"])))
    print(f"  Experiment {e}: {np.exp(lb[-1]):.2f} [{np.exp(lh[-1]):.2f}]")
for lab, v in (("Bliss", np.array(lb)), ("HSA", np.array(lh))):
    n = len(v); m, s = v.mean(), v.std(ddof=1); tq = stats.t.ppf(0.975, n - 1); p = 2 * stats.t.sf(abs(m / (s / np.sqrt(n))), n - 1)
    print(f"Across {n} experiments, observed / {lab} expected: {np.exp(m):.2f} (95% CI {np.exp(m - tq * s / np.sqrt(n)):.2f} to "
          f"{np.exp(m + tq * s / np.sqrt(n)):.2f}), one-sample t-test on the log ratios p = {p:.2g}")
