"""Authored raw-label logistic MLE calibrator and unique-score PAV; not the complete Platt algorithm or sklearn."""
from fractions import Fraction as F
from bisect import bisect_right
import math,json
S=[-3.,-2.,-1.,0.,1.,2.,3.,4.];Y=[0,0,1,0,1,0,1,1]
TEST=[(-2.5,0),(-.5,0),(.5,1),(2.5,1)]
def sigmoid(z):return 1/(1+math.exp(-z)) if z>=0 else math.exp(z)/(1+math.exp(z))
def nll(a,b):return sum(max(a*s+b,0)-y*(a*s+b)+math.log1p(math.exp(-abs(a*s+b))) for s,y in zip(S,Y))/len(S)
def fit_sigmoid():
    a,b=0.,0.;history=[nll(a,b)]
    for iteration in range(60):
        ps=[sigmoid(a*s+b) for s in S];ws=[p*(1-p) for p in ps]
        ga=sum(s*(p-y) for s,p,y in zip(S,ps,Y))/len(S);gb=sum(p-y for p,y in zip(ps,Y))/len(S)
        if max(abs(ga),abs(gb))<1e-10:break
        haa=sum(s*s*w for s,w in zip(S,ws))/len(S);hab=sum(s*w for s,w in zip(S,ws))/len(S);hbb=sum(ws)/len(S);det=haa*hbb-hab*hab
        if det<=0:raise ValueError('singular Hessian')
        da=(hbb*ga-hab*gb)/det;db=(-hab*ga+haa*gb)/det;step=1.
        while nll(a-step*da,b-step*db)>history[-1]+1e-14:
            step/=2
            if step<1e-10:raise ValueError('line search stalled')
        a-=step*da;b-=step*db;history.append(nll(a,b))
    assert a>0 and abs(b)>1e-3 and history[-1]<history[0]
    return a,b,history
def pav(scores,labels):
    if len(scores)!=len(labels) or len(set(scores))!=len(scores) or scores!=sorted(scores):raise ValueError('unique sorted scores required in this exercise')
    blocks=[]
    for i,y in enumerate(labels):
        blocks.append({'start':i,'end':i,'sum':F(y),'count':1})
        while len(blocks)>1 and blocks[-2]['sum']/blocks[-2]['count']>blocks[-1]['sum']/blocks[-1]['count']:
            r=blocks.pop();l=blocks.pop();blocks.append({'start':l['start'],'end':r['end'],'sum':l['sum']+r['sum'],'count':l['count']+r['count']})
    fitted=[F(0)]*len(scores)
    for b in blocks:
        for i in range(b['start'],b['end']+1):fitted[i]=b['sum']/b['count']
    assert all(a<=b for a,b in zip(fitted,fitted[1:]))
    return fitted,blocks
def demo():
    a,b,history=fit_sigmoid();fit,blocks=pav(S,Y);rows=[]
    for name,predict in [('raw_sigmoid',lambda s:sigmoid(s)),('fitted_sigmoid',lambda s:sigmoid(a*s+b)),('authored_left_step_isotonic',lambda s:float(fit[max(0,min(len(S)-1,bisect_right(S,s)-1))]))]:
        vp=[predict(s) for s in S];tp=[predict(s) for s,y in TEST]
        rows.append({'method':name,'validation_probabilities':vp,'heldout_probabilities':tp,'validation_brier':sum((p-y)**2 for p,y in zip(vp,Y))/len(Y),'heldout_brier':sum((p-y)**2 for p,(s,y) in zip(tp,TEST))/len(TEST)})
    return {'validation_scores':S,'validation_labels':Y,'heldout_rows':[{'score':s,'label':y} for s,y in TEST],'fitted_sigmoid_a':a,'fitted_sigmoid_b':b,'validation_logistic_nll_history':history,'isotonic_fitted_values':[str(v) for v in fit],'isotonic_blocks':[{**z,'sum':str(z['sum']),'mean':str(z['sum']/z['count'])} for z in blocks],'new_score_rule':'left step using nearest fitted score at or below query, clipped to endpoints; not sklearn interpolation','rows':rows,'base_classifier_trained':False,'heldout_used_for_fit':False,'platt_target_smoothing_or_multiclass_or_sklearn_run':False}
if __name__=='__main__':print(json.dumps(demo(),ensure_ascii=False,indent=2))