SPB Git

spb/wp29_uqo

Public
2.0 KB · 57 lines python
1"""F2-smoke-01: OLS slope of y on x with firm-clustered SE. Run from workspace root."""2import json, os3import numpy as np, pandas as pd, statsmodels.api as sm4from scipy import stats56df_raw = pd.read_csv("data/panel.csv")7n_input = len(df_raw)8issues = []9df = df_raw.copy()10na = df[["firm_id", "x", "y"]].isna().any(axis=1)11if na.any():12    issues.append({"issue": f"{int(na.sum())} rows with missing values", "handling": "dropped"})13    df = df[~na]14dup = df.duplicated()15if dup.any():16    issues.append({"issue": f"{int(dup.sum())} exact duplicate rows", "handling": "kept (cannot distinguish from genuine repeats); flagged"})17nonfinite = ~np.isfinite(df[["x", "y"]]).all(axis=1)18if nonfinite.any():19    issues.append({"issue": f"{int(nonfinite.sum())} non-finite values", "handling": "dropped"})20    df = df[~nonfinite]2122X = sm.add_constant(df["x"].to_numpy())23y = df["y"].to_numpy()24g = df["firm_id"].to_numpy()25ols = sm.OLS(y, X).fit()26cl = sm.OLS(y, X).fit(cov_type="cluster", cov_kwds={"groups": g, "use_correction": True})27G = int(pd.Series(g).nunique())28beta, alpha = float(ols.params[1]), float(ols.params[0])29se_cl = float(cl.bse[1])30tcrit = float(stats.t.ppf(0.975, G - 1))  # t with G-1 df for clustered inference3132res = {33    "task_id": "F2-smoke-01",34    "beta_hat": beta,35    "alpha_hat": alpha,36    "se_ols": float(ols.bse[1]),37    "se_cluster": se_cl,38    "beta_ci_low": beta - tcrit * se_cl,39    "beta_ci_high": beta + tcrit * se_cl,40    "n_analysed": int(len(df)),41    "n_clusters": G,42    "feasibility": "complete",43    "feasibility_reason": None,44    "data_issues": issues,45    "n_input_rows": {"panel.csv": int(n_input)},46    "random_seed": None,47    # supplementary diagnostics48    "se_cluster_ratio_to_ols": se_cl / float(ols.bse[1]),49    "t_crit_df": G - 1,50    "t_crit": tcrit,51    "cluster_sizes_min_max": [int(pd.Series(g).value_counts().min()), int(pd.Series(g).value_counts().max())],52}53os.makedirs("out", exist_ok=True)54with open("out/results.json", "w") as f:55    json.dump(res, f, indent=2)56print(json.dumps(res, indent=2))57
Files · main
Menu

Theme