# Calcium imaging responses, WT vs KO: neurons as n, mice as n, and the effect of selecting "responsive" neurons on the same
# trials that measure them. Python 3.13, NumPy, pandas, SciPy.
import numpy as np, pandas as pd
from scipy import stats

d = pd.read_csv("responses.csv")                     # one row per trial: mouse, genotype, neuron, trial, dff
nm = d.groupby(["mouse", "genotype", "neuron"]).dff.mean().reset_index()      # one response per neuron
ko = nm.genotype == "KO"
t1 = stats.ttest_ind(nm.dff[ko], nm.dff[~ko], equal_var=False)
print(f"Neurons as n ({ko.sum()} vs {(~ko).sum()} neurons): KO - WT {nm.dff[ko].mean() - nm.dff[~ko].mean():.3f}, Welch's t-test p = {t1.pvalue:.2g}")
mm = nm.groupby(["mouse", "genotype"]).dff.mean().reset_index(); km = mm.genotype == "KO"
t2 = stats.ttest_ind(mm.dff[km], mm.dff[~km], equal_var=False); c2 = t2.confidence_interval()
print(f"Mice as n (5 vs 5 mice, mean of each mouse's neurons): KO - WT {mm.dff[km].mean() - mm.dff[~km].mean():.3f} "
      f"(95% CI {c2.low:.3f} to {c2.high:.3f}), p = {t2.pvalue:.2g}")
# "Responsive" neurons: one-sample t-test over trials, p < 0.05 and mean above 0
odd, even = d[d.trial % 2 == 1], d[d.trial % 2 == 0]
def respond(x):
    g = x.groupby(["mouse", "genotype", "neuron"]).dff
    return pd.DataFrame({"mean": g.mean(), "p": g.apply(lambda v: stats.ttest_1samp(v, 0).pvalue)}).reset_index()
allt = respond(d); sel = (allt.p < 0.05) & (allt["mean"] > 0)
print(f"Responsive on all 20 trials: {sel.sum()} of {len(allt)} neurons ({sel.mean():.1%}); their mean ΔF/F on the same trials {allt['mean'][sel].mean():.3f}; "
      f"mean ΔF/F of all neurons {allt['mean'].mean():.3f} (true average 0.05)")
so = respond(odd); se = respond(even); s_odd = (so.p < 0.05) & (so["mean"] > 0)
print(f"Selected on odd trials ({s_odd.sum()} neurons): mean ΔF/F on odd trials {so['mean'][s_odd].mean():.3f}, on held-out even trials {se['mean'][s_odd].mean():.3f}")
pm = allt.assign(resp=sel).groupby(["mouse", "genotype"]).resp.mean().reset_index(); kp = pm.genotype == "KO"
t3 = stats.ttest_ind(pm.resp[kp], pm.resp[~kp], equal_var=False)
print(f"Share of responsive neurons per mouse: WT {pm.resp[~kp].mean():.1%}, KO {pm.resp[kp].mean():.1%}; mice as n, Welch's t-test p = {t3.pvalue:.2g}")
