import json
from pathlib import Path
import numpy as np
import pandas as pd
p=pd.read_csv('/input/data/analysis_panel.csv')
p=p.dropna(subset=['y_opioid_rate','law_direct','synthetic_dominant'])
x=np.column_stack([np.ones(len(p)),p.law_direct,p.synthetic_dominant,p.law_direct*p.synthetic_dominant,pd.get_dummies(p.state,drop_first=True,dtype=float),pd.get_dummies(p.year,drop_first=True,dtype=float)])
y=p.y_opioid_rate.to_numpy()
b=np.linalg.lstsq(x,y,rcond=None)[0]
res=y-x@b; bread=np.linalg.pinv(x.T@x); meat=np.zeros((x.shape[1],x.shape[1]))
for state in sorted(p.state.unique()):
 mask=(p.state==state).to_numpy();score=x[mask].T@res[mask];meat+=np.outer(score,score)
n,k=x.shape;g=p.state.nunique();cov=(g/(g-1))*((n-1)/(n-k))*bread@meat@bread
se=np.sqrt(np.diag(cov));z=1.959963984540054
out={'n':n,'n_states':int(g),'rank':int(np.linalg.matrix_rank(x)),'law':float(b[1]),'law_ci95':[float(b[1]-z*se[1]),float(b[1]+z*se[1])],'interaction':float(b[3]),'dominant':float(b[1]+b[3])}
print(json.dumps(out,indent=2))
