"""Small exact-rational least squares in a polynomial basis; no numerical-library fit."""
from fractions import Fraction as F
import json
def solve(matrix,rhs):
    n=len(rhs);a=[[F(x) for x in row]+[F(y)] for row,y in zip(matrix,rhs)]
    for c in range(n):
        pivot=next((i for i in range(c,n) if a[i][c]),None)
        if pivot is None:raise ValueError('singular normal matrix in this limited solver')
        a[c],a[pivot]=a[pivot],a[c];scale=a[c][c];a[c]=[x/scale for x in a[c]]
        for i in range(n):
            if i!=c:
                t=a[i][c];a[i]=[x-t*y for x,y in zip(a[i],a[c])]
    return [a[i][-1] for i in range(n)]
def fit(xs,ys,degree):
    design=[[F(x)**j for j in range(degree+1)] for x in xs]
    p=degree+1;normal=[[sum(row[j]*row[k] for row in design) for k in range(p)] for j in range(p)]
    rhs=[sum(row[j]*F(y) for row,y in zip(design,ys)) for j in range(p)]
    weights=solve(normal,rhs);pred=[sum(v*b for v,b in zip(row,weights)) for row in design]
    residual=[F(y)-p for y,p in zip(ys,pred)]
    assert all(sum(row[j]*e for row,e in zip(design,residual))==0 for j in range(p))
    return {'degree':degree,'design':[[str(x) for x in row] for row in design],'coefficients':[str(x) for x in weights],'predictions':[str(x) for x in pred],'residual_y_minus_pred':[str(x) for x in residual],'training_SSE':str(sum(x*x for x in residual)),'x3_extrapolation':str(sum(b*3**j for j,b in enumerate(weights)))}
def demo():
    xs=[-2,-1,0,1,2];ys=[5,2,1,2,5];linear=fit(xs,ys,1);quadratic=fit(xs,ys,2)
    assert linear['coefficients']==['3','0'] and linear['training_SSE']=='14'
    assert quadratic['coefficients']==['1','0','1'] and quadratic['training_SSE']=='0'
    assert linear['x3_extrapolation']=='3' and quadratic['x3_extrapolation']=='10'
    # Linearity in coefficients: phi(x)*(u+v) == phi(x)*u + phi(x)*v.
    u=[F(1),F(2),F(3)];v=[F(2),F(-1),F(1)];phi=[F(1),F(2),F(4)]
    assert sum(p*(a+b) for p,a,b in zip(phi,u,v))==sum(p*a for p,a in zip(phi,u))+sum(p*b for p,b in zip(phi,v))
    return {'authored_x':xs,'authored_y':ys,'linear':linear,'quadratic':quadratic,'coefficient_linearity_checked':True,'outside_x3_observed_label':None,'holdout_or_generalization_test':False,'floating_point_solver_recommendation':False}
if __name__=='__main__':print(json.dumps(demo(),ensure_ascii=False,indent=2))
