1"""Calibration of the 2026 base point (US) from data/processed/calibration_targets.csv.23Units: GDP in 2026 = 1 (per efficiency unit of labor/land/resource, gf = 1); workers H_L + H_H = 1.4The base point is a balanced-growth path of the frozen-AI economy (alpha fixed at its 2026 value). CES weights5are set so that the base-point cost shares equal the targets (share_j = nu_j x_j^rho / sum, hence6nu_j ∝ share_j x_j^{-rho}); the aggregate scale Z_0 then normalizes final output. Every target has a7`source` and a `status` ("data" if read from the data file, "assumed" if a documented assumption,8"scenario" if a forward-looking scenario input). The result is checked by re-solving the BGP (steady.py).9"""10from __future__ import annotations1112import csv13import json14import os15from dataclasses import dataclass, field, asdict16from typing import Dict1718import numpy as np1920from .ces import ces, rho_of_sigma21from .params import Params, TaskSpec22from .steady import bgp_interest_rate, solve_bgp, target_rentals2324ROOT = os.path.abspath(os.path.join(os.path.dirname(__file__), "..", ".."))25TARGETS_CSV = os.path.join(ROOT, "data", "processed", "calibration_targets.csv")262728@dataclass29class Target:30 value: float31 source: str32 status: str = "assumed" # data | assumed | scenario333435def _default_targets() -> Dict[str, Target]:36 """Documented assumptions used when a data target is not retrieved. Overwritten by the data file."""37 A = "assumed (see docs/assumptions_register.md)"38 L = "literature value (see docs/assumptions_register.md)"39 return {40 "labor_share": Target(0.58, A),41 "housing_va_share": Target(0.12, A + "; BEA housing services ~12% of GDP, not retrieved as a series"),42 "KN_over_GDP": Target(1.22, A), # overwritten by data; K_N is residual anyway43 "KH_over_GDP": Target(1.17, A),44 "KE_over_GDP": Target(0.30, A),45 "VQ_over_GDP": Target(0.06, A), # market value of compute capital / GDP (derived, see CSV_DERIVED)46 "energy_share": Target(0.04, A), # energy sector gross output / GDP47 "datacenter_energy_share": Target(0.044, A), # share of electricity used by data centers48 "land_rent_share": Target(0.03, A), # production land rent / GDP49 "resource_share_energy": Target(0.30, A), # eta (implied in calibration)50 "housing_land_share": Target(0.35, A), # nu_H (implied in calibration)51 "delta_N": Target(0.089, A),52 "delta_H": Target(0.023, A),53 "delta_E": Target(0.05, A),54 "delta_Q": Target(0.25, A),55 "r_star": Target(0.035, A), # real return on wealth in the frozen-AI BGP56 "g_Z": Target(0.015, A), # trend growth per efficiency unit57 "gamma_u": Target(2.0, L),58 "share_L_workers": Target(0.60, A), # share of workers with skill L59 "skill_premium": Target(1.7, A), # w_H / w_L per worker60 "task_wage_shares_L": Target(0.40, A), # share of L wage bill in task 1 (rest in task 2)61 "task_wage_shares_H": Target(0.50, A), # share of H wage bill in task 3 (rest in task 4)62 "sigma_B": Target(0.5, L),63 "sigma_N": Target(0.7, L),64 "sigma_T": Target(0.6, L),65 "phi_N": Target(3.0, L), "phi_Q": Target(1.0, L), "phi_E": Target(3.0, L), "phi_H": Target(3.0, L),66 "g_pQ": Target(0.15, "scenario"), "g_eps": Target(0.10, "scenario"), "g_alpha": Target(0.10, "scenario"),67 "pQ_floor_ratio": Target(1e-2, "scenario"), "eps_floor_ratio": Target(1e-2, "scenario"),68 }697071# calibration_targets.csv target_id -> target name (value used directly)72CSV_MAP = {73 "T01_labor_share_PWT_US_latest": "labor_share",74 "T06_net_stock_private_nonres_fixed_assets_GDP_latest": "KN_over_GDP",75 "T06_depreciation_rate_private_nonres_latest": "delta_N",76 "T10_datacenter_share_of_US_electricity_2023": ("datacenter_energy_share", 0.01), # percent -> ratio77}7879# derived data targets: name -> (function of the CSV rows, description)80def _derived(rows):81 out = {}82 g = lambda k: float(rows[k]["value"])83 try:84 tot, nonres = g("T06_net_stock_private_fixed_assets_GDP_latest"), g("T06_net_stock_private_nonres_fixed_assets_GDP_latest")85 out["KH_over_GDP"] = Target(tot - nonres, "BEA fixed assets via FRED: total minus nonresidential private net stock / GDP (2024)", "data")86 dtot, dnr = g("T06_depreciation_rate_private_total_latest"), g("T06_depreciation_rate_private_nonres_latest")87 out["delta_H"] = Target((dtot * tot - dnr * nonres) / (tot - nonres), "implied residential depreciation = (CFC_total - CFC_nonres)/K_res (BEA via FRED, 2024)", "data")88 ipe = g("T05_info_processing_equipment_GDP_latest")89 out["VQ_over_GDP"] = Target(ipe / (0.25 + 0.015 + 0.15), "info-processing-equipment investment/GDP (FRED Y034RC1Q027SBEA, 2026Q2) capitalized at delta_Q + g_Z + g_pQ", "data-derived")90 except KeyError:91 pass92 return out939495def load_targets(path=TARGETS_CSV) -> Dict[str, Target]:96 tg = _default_targets()97 if os.path.exists(path):98 with open(path, newline="") as f:99 rows = {r["target_id"]: r for r in csv.DictReader(f)}100 for tid, spec in CSV_MAP.items():101 row = rows.get(tid)102 if row and row.get("status") == "RETRIEVED" and row.get("value") not in (None, ""):103 key, scale = (spec, 1.0) if isinstance(spec, str) else spec104 tg[key] = Target(float(row["value"]) * scale, f"{row.get('source')} {row.get('series_id')} ({row.get('year_or_period')})", "data")105 tg.update(_derived(rows))106 return tg107108109def calibrate(tg: Dict[str, Target] | None = None, verbose=False, overrides: Dict[str, float] | None = None):110 """overrides: target values to replace (e.g. scenario inputs g_pQ, sigma_N, r_star); the base point is111 re-derived so that 2026 remains consistent with the data under the scenario's anticipated price path."""112 tg = dict(load_targets() if tg is None else tg)113 for k, val in (overrides or {}).items():114 tg[k] = Target(val, f"override ({tg[k].source})" if k in tg else "override", tg[k].status if k in tg else "scenario")115 v = {k: t.value for k, t in tg.items()}116 g = v["g_Z"]117 gamma_u = v["gamma_u"]118 r = v["r_star"]119 beta = (1 + g) ** gamma_u / (1 + r)120 GDP = 1.0121 hv = v["housing_va_share"] * GDP122 Y = GDP - hv123 delta = {"N": v["delta_N"], "Q": v["delta_Q"], "E": v["delta_E"], "H": v["delta_H"]}124 phi = {"N": v["phi_N"], "Q": v["phi_Q"], "E": v["phi_E"], "H": v["phi_H"]}125 K = {"N": v["KN_over_GDP"] * GDP, "H": v["KH_over_GDP"] * GDP, "E": v["KE_over_GDP"] * GDP}126 VQ = v["VQ_over_GDP"] * GDP127 # --- rentals at the BGP (per unit of price); compute: user cost with the anticipated price decline128 # R_Q,0 = p_0 [ (1+r) e^{g_p} - (1 - delta_Q) ] (capital bought at t=-1 at price p_0 e^{g_p})129 Rst = target_rentals(Params(delta=delta, phi=phi, g_Z=g, beta=beta, gamma_u=gamma_u, pQ0=1.0), r, 1.0)130 Rst["Q"] = (1 + r) * np.exp(v["g_pQ"]) - (1 - delta["Q"])131 # --- energy block: E_0 = 1 ; P_E E = energy_share132 energy_val = v["energy_share"] * GDP133 E0 = 1.0134 PE = energy_val / E0135 RE = 1.0136 eta_implied = 1 - Rst["E"] * K["E"] / energy_val137 if verbose:138 print(f" resource share in energy implied by K_E target: {eta_implied:.3f} (assumed {v['resource_share_energy']:.3f}) -> using implied")139 eta = float(np.clip(eta_implied, 0.05, 0.9))140 ZE = E0 / (K["E"] ** (1 - eta) * RE ** eta)141 # --- labor: workers, wages142 HL = v["share_L_workers"]143 HH = 1 - HL144 wage_bill = v["labor_share"] * GDP145 wL = wage_bill / (HL + HH * v["skill_premium"])146 wH = wL * v["skill_premium"]147 # --- tasks: wage shares across tasks within skill; task 1 also uses compute148 sL1 = v["task_wage_shares_L"]149 sH3 = v["task_wage_shares_H"]150 LB = {0: sL1 * HL, 1: (1 - sL1) * HL, 2: sH3 * HH, 3: (1 - sH3) * HH} # labor in each task151 # --- compute block: K_Q = 1 by normalization of compute units; p^Q_0 = VQ; the compute Euler gives R_Q;152 # the services price A = R_Q + eps_0 P_E; the task-1 tie A = w_L/alpha_1 gives alpha_1(2026); the153 # compute share of task 1 is then implied (reported, not targeted).154 KQ = 1.0155 pQ0 = VQ / KQ156 eps0 = v["datacenter_energy_share"] * E0 / KQ157 RQ = Rst["Q"] * pQ0158 A = RQ + eps0 * PE159 alpha10 = wL / A160 Q_eff = KQ / alpha10 # compute in task 1, worker equivalents161 c1 = Q_eff / (LB[0] + Q_eff)162 if verbose:163 print(f" implied compute share of task-1 input in 2026: {c1:.3f}; alpha_1(2026) = {alpha10:.1f}")164 EY = E0 - eps0 * KQ165 # --- structures: K_N is residual so that factor payments add up to GDP = 1 exactly166 land_val = v["land_rent_share"] * GDP167 energy_in_Y = PE * (E0 - eps0 * KQ)168 RN = Rst["N"]169 valN = Y - wage_bill - A * KQ - land_val - energy_in_Y170 K["N"] = valN / RN171 if verbose:172 print(f" K_N/GDP implied by income shares: {K['N']:.3f} (target {v['KN_over_GDP']:.3f})")173 # --- consumption / housing174 inv = sum((delta[k] + g) * K[k] + 0.5 * phi[k] * g ** 2 * K[k] for k in K) + pQ0 * ((delta["Q"] + g) * KQ + 0.5 * phi["Q"] * g ** 2 * KQ)175 C = Y - inv176 theta_H = hv / (hv + C)177 X = np.array([LB[0] + Q_eff, LB[1], LB[2], LB[3]])178 Fval = np.array([wL, wL, wH, wH]) # marginal values of tasks179 task_val = Fval * X # = wages + compute income in task 1180 rhoB = rho_of_sigma(v["sigma_B"])181 nuB = task_val * X ** (-rhoB)182 nuB = nuB / nuB.sum()183 XB = ces(nuB, X, rhoB)184 PB = task_val.sum() / XB # price of the composite185 # --- V nest: structures vs tasks186 rhoN = rho_of_sigma(v["sigma_N"])187 valN, valB = RN * K["N"], PB * XB188 wN = np.array([valN, valB]) * np.array([K["N"], XB]) ** (-rhoN)189 wN = wN / wN.sum()190 nu_N = float(wN[0])191 V = ces([nu_N, 1 - nu_N], [K["N"], XB], rhoN)192 PV = (valN + valB) / V193 # --- top nest: land, energy, V ; land rent share target194 rhoT = rho_of_sigma(v["sigma_T"])195 T_land = 1.0196 valE = PE * EY197 # Y (gross output of final goods) must equal valT + valE + valV + accessory(=0 at base)198 valV = valN + valB199 Y_implied = land_val + valE + valV200 wT = np.array([land_val, valE, valV]) * np.array([T_land, EY, V]) ** (-rhoT)201 wT = wT / wT.sum()202 nu_T, nu_E = float(wT[0]), float(wT[1])203 Gval = ces(wT, [T_land, EY, V], rhoT)204 Z0 = Y_implied / Gval205 # --- housing: land share implied by the structures' rental and the housing value-added target206 nu_H_implied = 1 - Rst["H"] * K["H"] / hv207 if verbose:208 print(f" land share in housing implied by K_H target: {nu_H_implied:.3f} (assumed {v['housing_land_share']:.3f}) -> using implied")209 nu_H = float(np.clip(nu_H_implied, 0.05, 0.9))210 TH = 1.0211 # --- tasks of the quantitative model: four bottleneck classes. Accessory tasks (linear add-on) are kept212 # for the Restrepo replication only: with an energy constraint a linear accessory block makes the213 # value of compute bounded away from zero and lets compute absorb all energy (docs/research_log.md).214 tasks = TaskSpec(names=["routine", "cognitive-L", "cognitive-H", "expert"],215 skill=[0, 0, 1, 1], nu=nuB.tolist(), gamma=[0.0, 0.0, 0.0, 0.0],216 alpha_inf=[1.0, 2.0, 4.0, 8.0])217 alpha0_mult = [alpha10 / 1.0, 50.0, 400.0, 3000.0] # task 1 anchored to 2026 data; others are scenario inputs218 p = Params(tasks=tasks, sigma_B=v["sigma_B"], H=[HL, HH], sigma_N=v["sigma_N"], nu_N=nu_N, sigma_T=v["sigma_T"],219 nu_T=nu_T, nu_E=nu_E, Z0=Z0, g_Z=g, ZE=ZE, eta=eta, RE=RE, T_land=T_land, TH_land=TH, nu_H=nu_H,220 theta_H=theta_H, delta=delta, phi=phi, pQ0=pQ0, g_pQ=v["g_pQ"], pQ_floor_ratio=v["pQ_floor_ratio"],221 eps0=eps0, g_eps=v["g_eps"], eps_floor_ratio=v["eps_floor_ratio"], beta=beta, gamma_u=gamma_u,222 alpha0_mult=alpha0_mult, g_alpha=v["g_alpha"], KN0=K["N"], KQ0=KQ, KE0=K["E"], KH0=K["H"])223 p.source = {k: f"{t.status}: {t.source}" for k, t in tg.items()}224 p.source["eta"] = "calibrated: implied by K_E/GDP and energy share"225 p.source["KN0"] = "calibrated: residual so that factor payments add up to GDP (K_N/GDP target not imposed)"226 p.source["nu_H"] = "calibrated: implied by K_H/GDP, delta_H, r* and the housing value-added share"227 p.source["theta_H"] = "calibrated: housing VA share and goods consumption at the base point"228 p.source["alpha_1_2026"] = "calibrated: value of compute capital, compute Euler and task-1 tie A = w_L/alpha_1"229 p.source["compute_in_task1_implied"] = f"{c1:.4f}"230 p.source["Y_implied"] = f"{Y_implied:.4f} (final goods gross output; GDP = 1 by normalization)"231 info = dict(Y=Y_implied, C=C, theta_H=theta_H, A=A, alpha10=alpha10, KQ=KQ, pQ0=pQ0, eps0=eps0, wL=wL, wH=wH,232 nuB=nuB.tolist(), nu_N=nu_N, nu_T=nu_T, nu_E=nu_E, Z0=Z0, eta=eta, beta=beta, r=r, RQ=RQ, PE=PE)233 return p, info234235236def verify_calibration(p, info, tol=1e-6, verbose=True):237 """Solve the frozen-AI BGP from the calibrated stocks: the stocks must be a fixed point."""238 b = solve_bgp(p.replace(alpha_freeze=True), verbose=verbose)239 dev = {k: abs(np.log(b["K"][k] / {"N": p.KN0, "Q": p.KQ0, "E": p.KE0, "H": p.KH0}[k])) for k in b["K"]}240 dev["pQ0"] = abs(np.log(b["pQ0"] / p.pQ0))241 per = b["period"]242 out = dict(max_dev=max(dev.values()), dev=dev, GDP=per["GDP"], labor_share=per["wage_bill"] / per["GDP"],243 compute_income_share=per["RQ"] * per["KQ"] / per["GDP"], energy_share=per["PE"] * per["E"] / per["GDP"],244 land_rent_share=per["RT"] * per["T_eff"] / per["GDP"], housing_share=per["RH"] * per["HS"] / per["GDP"],245 A=per["A"], util=per["util"], wL=per["w"][0], wH=per["w"][1])246 if verbose:247 print(" calibration check:", {k: (round(float(x), 5) if not isinstance(x, dict) else x) for k, x in out.items()})248 return out, b249250251def save_parameters(p, info, path):252 d = p.to_dict()253 d["_calibration_info"] = {k: (v if not hasattr(v, "tolist") else v.tolist()) for k, v in info.items()}254 with open(path, "w") as f:255 json.dump(d, f, indent=2, default=float)256