# How often each analysis gives p < 0.05: 5,000 simulated studies of 4 experiments, with no knockdown effect and with
# a knockdown that cuts migration by 30%, and the 30% cut with 8 experiments. Python 3.13, NumPy, SciPy. Seed 20261112.
import numpy as np, warnings
from scipy import stats
warnings.simplefilter("ignore", RuntimeWarning)                    # SciPy warns that the fixed 100% column has no variance
from model import experiment, INSERTS, FIELDS
rng = np.random.default_rng(20261112); R = 5_000

def one(kd_ratio, E):
    c = np.array([[r[3] for r in experiment(rng, kd_ratio=kd_ratio)] for _ in range(E)], float).reshape(E, 2, INSERTS, FIELDS)
    ins = c.mean(axis=3); ex = ins.mean(axis=2)                          # insert means; experiment means [experiment, condition]
    pa = stats.ttest_ind(c[:, 1].ravel(), c[:, 0].ravel(), equal_var=False).pvalue
    pb = stats.ttest_ind(ins[:, 1].ravel(), ins[:, 0].ravel(), equal_var=False).pvalue
    pc = stats.ttest_ind(100 * ex[:, 1] / ex[:, 0], np.full(E, 100.0)).pvalue
    pd_ = stats.ttest_1samp(100 * ex[:, 1] / ex[:, 0], 100).pvalue
    pe = stats.ttest_1samp(np.log(ex[:, 1] / ex[:, 0]), 0).pvalue
    return np.array([pa, pb, pc, pd_, pe]) < 0.05
for lab, kr, E in (("Knockdown has no effect, 4 experiments", 1.0, 4), ("Knockdown cuts migration by 30%, 4 experiments", 0.7, 4),
                   ("Knockdown cuts migration by 30%, 8 experiments", 0.7, 8)):
    r = np.array([one(kr, E) for _ in range(R)]).mean(axis=0)
    print(f"{lab}: p<0.05 with A fields as n {r[0]:.1%}; B inserts as n {r[1]:.1%}; C percent of control, control fixed at 100% {r[2]:.1%}; "
          f"D percent of control, one-sample vs 100 {r[3]:.1%}; E log ratio per experiment {r[4]:.1%}")