import math, itertools
from statistics import NormalDist
# t quantile via scipy-free: use known t_{3,0.975}=3.182446305, t_{4}=2.776445105
T={3:3.182446305284263,4:2.7764451051977987,2:4.302652729911275}
def meta(tr):
    y=[math.log((a/n)/(c/m)) for a,n,c,m in tr]; v=[1/a-1/n+1/c-1/m for a,n,c,m in tr]
    w=[1/x for x in v]; k=len(y); mu=sum(wi*yi for wi,yi in zip(w,y))/sum(w)
    Q=sum(wi*(yi-mu)**2 for wi,yi in zip(w,y)); t2=max(0,(Q-k+1)/(sum(w)-sum(x*x for x in w)/sum(w)))
    ws=[1/(x+t2) for x in v]; m2=sum(wi*yi for wi,yi in zip(ws,y))/sum(ws)
    s=max(1,sum(wi*(yi-m2)**2 for wi,yi in zip(ws,y))/(k-1)); se=math.sqrt(s/sum(ws))
    return math.exp(m2),math.exp(m2-T[k-1]*se),math.exp(m2+T[k-1]*se),t2
base=[(34,2157,27,1806),(10,817,9,774),(7,708,2,722)]
co=(4,374,5,356)
print("prev",meta(base)); print("primary",meta(base+[co]))
print("n",sum(a[1]+a[3] for a in base+[co]))
# together unions
cnt=0;inc=0;rr=[];lo=[];hi=[]
for ut in range(79,min(679,79+21)+1):
  for uc in range(95,min(679,95+24)+1):
    r=meta(base+[co,(ut,679,uc,679)]);cnt+=1;inc+= r[1]<=1<=r[2];rr.append(r[0]);lo.append(r[1]);hi.append(r[2])
print("together",cnt,inc,min(rr),max(rr),min(lo),max(hi))
# denominator
rs=[]
for et in range(0,37):
  for ec in range(0,43):
    rs.append(((4+et)/410)/((5+ec)/398))
print("denom",len(rs),min(rs),max(rs))
for et in range(37):
  if ((4+et)/410)/(5/398)>1: print("reverse at",et);break
print("registry-own-denominator gaps",407-374,391-356)