import json,csv,math,statistics,sys
from collections import Counter
from pathlib import Path
p=Path(sys.argv[1]); a=list(json.loads((p/'results/associations.json').read_text()).values())
def bh(xs):
 n=len(xs);o=sorted(range(n),key=xs.__getitem__);q=[0.]*n;v=1.
 for j in range(n-1,-1,-1):
  i=o[j];v=min(v,xs[i]*n/(j+1));q[i]=v
 return q
for field,out in [('replication','q'),('published_difference_p','published_difference_q'),('published_difference_p_t','published_difference_q_t'),('heterogeneity_p','heterogeneity_q'),('heterogeneity_p_t','heterogeneity_q_t')]:
 xs=[r[field]['p'] if field=='replication' else r[field] for r in a];q=bh(xs)
 assert all(abs(v-r[out])<1e-11 for r,v in zip(a,q)),out
rep=[r['replication']['b']*r['published_b'][0]>0 and q<.05 for r,q in zip(a,bh([r['replication']['p'] for r in a]))]
repro=[r['published_b'][1]<=r['original']['b']<=r['published_b'][2] for r in a]
res={'n':len(a),'replicated':sum(rep),'reproduced':sum(repro),'nominal_power80':sum(r['power']>=.8 for r in a),'replicated_power80':sum(v and r['power']>=.8 for r,v in zip(a,rep)),'bonferroni_power80':sum(r['power_bonferroni']>=.8 for r in a),'replicated_bonferroni_power80':sum(v and r['power_bonferroni']>=.8 for r,v in zip(a,rep)),'inside_published_CI':sum(r['published_b'][1]<=r['replication']['b']<=r['published_b'][2] for r in a),'same_sign_nominal':sum(r['replication']['b']*r['published_b'][0]>0 and r['replication']['p']<.05 for r in a),'replicated_reproduced':sum(v and w for v,w in zip(rep,repro))}
for name in ['published_difference_q','published_difference_q_t','heterogeneity_q','heterogeneity_q_t']:res[name+'_under05']=sum(r[name]<.05 for r in a)
for name,num,den in [('replication_to_published','replication',None),('original_to_published','original',None),('replication_to_harmonized','replication','harmonized')]:res['median_'+name]=statistics.median(r[num]['b']/(r[den]['b'] if den else r['published_b'][0]) for r in a)
rows=list(csv.DictReader((p/'data/departure_coding.csv').open()));n=len(rows);obs=sum(r['coder_1']==r['coder_2'] for r in rows)/n;c1=Counter(r['coder_1'] for r in rows);c2=Counter(r['coder_2'] for r in rows);chance=sum(c1[k]*c2[k] for k in set(c1)|set(c2))/n**2
res['coding']={'n':n,'agree':round(obs*n),'kappa':(obs-chance)/(1-chance)}
res['headline_departures']=dict(Counter(d['evidence'] for r in a for d in r['departures'] if d['affects_headline']))
res['scope']='Independent arithmetic from supplied aggregate fit outputs, not fresh-data reproduction; no participant records used.'
print(json.dumps(res,indent=2))
