1"""F5-conf-01: Walrasian equilibrium of a 2x2 Cobb-Douglas exchange economy.2Run from workspace root: ~/wp29/venv/bin/python out/analysis.py3Reads only data/; writes out/results.json deterministically (no randomness).4"""5import json6from pathlib import Path7import numpy as np8import pandas as pd9from scipy.optimize import brentq1011ROOT = Path(__file__).resolve().parent.parent12end = pd.read_csv(ROOT / "data" / "endowments.csv")13par = pd.read_csv(ROOT / "data" / "parameters.csv")14P = dict(zip(par["parameter"], par["value"].astype(float)))1516issues = []17# --- validation -----------------------------------------------------------18assert set(end["agent"]) == {"A", "B"} and not end["agent"].duplicated().any()19assert end[["good1", "good2"]].notna().all().all()20assert (end[["good1", "good2"]] > 0).all().all(), "endowments must be strictly positive"21aA, aB = P["alpha_A"], P["alpha_B"]22assert 0 < aA < 1 and 0 < aB < 123e = end.set_index("agent")[["good1", "good2"]].astype(float)24alpha = {"A": aA, "B": aB}25agents = ["A", "B"]2627W1 = float(e["good1"].sum())28W2 = float(e["good2"].sum())29for g, W, key in [("good1", W1, "aggregate_endowment_good1"), ("good2", W2, "aggregate_endowment_good2")]:30 rec = P[key]31 if not np.isclose(rec, W):32 issues.append({33 "issue": f"parameters.csv {key} = {rec:g} disagrees with the sum of individual endowments in endowments.csv ({W:g}).",34 "handling": f"Used {W:g} (sum of primary individual endowments) so that aggregate supply is consistent with "35 f"agents' budget constraints; summary-sheet value {rec:g} discarded (reported only as a sensitivity in report.md).",36 })3738# --- demand / excess demand ------------------------------------------------39def demand(i, p):40 m = p * e.loc[i, "good1"] + e.loc[i, "good2"]41 return alpha[i] * m / p, (1 - alpha[i]) * m4243def z(p, W1=W1, W2=W2):44 x1 = sum(demand(i, p)[0] for i in agents)45 x2 = sum(demand(i, p)[1] for i in agents)46 return x1 - W1, x2 - W24748# --- closed form: sum_i a_i (p e1_i + e2_i) = p W1 => p = sum a_i e2_i / sum (1-a_i) e1_i49p_cf = sum(alpha[i] * e.loc[i, "good2"] for i in agents) / sum((1 - alpha[i]) * e.loc[i, "good1"] for i in agents)5051# --- independent numerical check: root of z1(p) via Brent on a wide bracket -52p_num = brentq(lambda p: z(p)[0], 1e-6, 1e6, xtol=1e-14, rtol=1e-15, maxiter=500)5354x = {i: demand(i, p_cf) for i in agents}55z1, z2 = z(p_cf)56walras = p_cf * z1 + z257# Walras' law holds identically: check on a grid of off-equilibrium prices too58grid = np.logspace(-2, 2, 41)59walras_grid_max = float(max(abs(p * z(p)[0] + z(p)[1]) for p in grid))6061# sensitivity: price that would clear good 1 if the summary-sheet aggregate were used62sens = {}63if not np.isclose(P["aggregate_endowment_good1"], W1):64 Ws = P["aggregate_endowment_good1"]65 ps = brentq(lambda p: z(p, W1=Ws)[0], 1e-6, 1e6, xtol=1e-14)66 zs1, zs2 = z(ps, W1=Ws)67 sens = {"sensitivity_price_if_summary_aggregate_good1": ps,68 "sensitivity_walras_residual_if_summary_aggregate_good1": ps * zs1 + zs2}6970res = {71 "task_id": "F5-conf-01",72 "price_ratio_p1_over_p2": p_cf,73 "x1_A": x["A"][0], "x2_A": x["A"][1], "x1_B": x["B"][0], "x2_B": x["B"][1],74 "excess_demand_good1": z1,75 "excess_demand_good2": z2,76 "walras_law_residual": walras,77 "walras_law_residual_max_abs_offeq_grid": walras_grid_max,78 "aggregate_endowment_good1": W1,79 "aggregate_endowment_good2": W2,80 "price_ratio_numerical": p_num,81 "abs_diff_closed_form_vs_numerical": abs(p_cf - p_num),82 "income_A": p_cf * e.loc["A", "good1"] + e.loc["A", "good2"],83 "income_B": p_cf * e.loc["B", "good1"] + e.loc["B", "good2"],84 **sens,85 "feasibility": "complete",86 "feasibility_reason": None,87 "data_issues": issues,88 "n_input_rows": {"endowments.csv": int(len(end)), "parameters.csv": int(len(par))},89 "n_analysed": len(agents),90 "random_seed": None,91}92res = {k: (float(v) if isinstance(v, (np.floating, np.integer)) and not isinstance(v, bool) else v) for k, v in res.items()}93(ROOT / "out" / "results.json").write_text(json.dumps(res, indent=2) + "\n")94print(f"p closed form = {p_cf:.12g}; p numerical = {p_num:.12g}; z1 = {z1:.3g}; Walras residual = {walras:.3g}")95