Tutorial: Godley-Lavoie Model PC — TIR Portfolio Choice¶

Model PC (Godley & Lavoie 2007, Chapter 4) is the canonical Stock-Flow Consistent benchmark: five sectors, two financial assets (money and bills), Tobin portfolio choice. It is the simplest model in which households actively allocate wealth between a zero-yield asset (money) and an interest-bearing asset (government bills).

EconIAC replaces the GL Tobin coefficients with TIR routing at rationality β:

$$w = \text{gibbs\_weights}([U_{\text{money}}, U_{\text{bills}}],\, \beta)$$ $$H_h = w_0 \cdot V \qquad B_h = w_1 \cdot V$$

where $V = H_h + B_h$ is household wealth and the utilities are $U_{\text{money}} = 0$, $U_{\text{bills}} = r$. This is a one-parameter family indexed by $\beta$:

β	Portfolio	Interpretation
0	50% money, 50% bills	Maximum entropy — indifferent
finite	Gibbs-weighted	Partially rational
→ ∞	0% money, 100% bills	Fully rational — classical Tobin

What this tutorial demonstrates:

Replicate GL Table 4.4 macro aggregates (Y, C, T, YD)
Time series: convergence to steady state
The portfolio share curve s(β) — the calibration landscape
β* calibration: recover implied rationality from national accounts data
The UK anomaly: s > 0.5 and what it means

Reference: Godley & Lavoie (2007) Monetary Economics, Chapter 4.
EconIAC implementation: Buckley (2026) TIR. doi:10.5281/zenodo.20237288

In [ ]:

Copied!





# !pip install econiac "jax[cpu]"

import jax.numpy as jnp
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import numpy as np

from econiac.core.ensemble import gibbs_weights
from econiac.economics.gl_pc import (
    PCParameters, ModelPC, calibrate_beta,
    portfolio_share_curve, GL_STEADY_STATE, GL_MONEY_SHARE,
)

plt.rcParams.update({"figure.dpi": 120, "font.size": 11})
# !pip install econiac "jax[cpu]"

import jax.numpy as jnp
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import numpy as np

from econiac.core.ensemble import gibbs_weights
from econiac.economics.gl_pc import (
    PCParameters, ModelPC, calibrate_beta,
    portfolio_share_curve, GL_STEADY_STATE, GL_MONEY_SHARE,
)

plt.rcParams.update({"figure.dpi": 120, "font.size": 11})

1. Calibrate β to GL Table 4.4¶

Godley & Lavoie's published steady state has $H_h = 21.62$, $V = 80$, so the money share is $s = H_h/V \approx 0.270$. We recover the β that the GL parameterisation implicitly assumes by minimising:

$$\mathcal{L}(\beta) = \bigl(\text{gibbs\_weights}([0, r],\,\beta)[0] - s_{\text{obs}}\bigr)^2$$

via JAX gradient descent.

In [ ]:

Copied!





beta_gl, losses_gl = calibrate_beta(
    GL_MONEY_SHARE, r_bar=0.025, n_steps=2000, lr=2.0, beta_init=10.0,
)
print(f"GL implied β*  = {beta_gl:.2f}")
print(f"Final loss     = {losses_gl[-1]:.2e}")
print(f"Target share   = {GL_MONEY_SHARE:.4f}")
w_pred = gibbs_weights(jnp.array([0.0, 0.025]), beta=beta_gl)
print(f"Predicted share = {float(w_pred[0]):.4f}  (error = {abs(float(w_pred[0]) - GL_MONEY_SHARE)*100:.4f}%)")
beta_gl, losses_gl = calibrate_beta(
    GL_MONEY_SHARE, r_bar=0.025, n_steps=2000, lr=2.0, beta_init=10.0,
)
print(f"GL implied β*  = {beta_gl:.2f}")
print(f"Final loss     = {losses_gl[-1]:.2e}")
print(f"Target share   = {GL_MONEY_SHARE:.4f}")
w_pred = gibbs_weights(jnp.array([0.0, 0.025]), beta=beta_gl)
print(f"Predicted share = {float(w_pred[0]):.4f}  (error = {abs(float(w_pred[0]) - GL_MONEY_SHARE)*100:.4f}%)")

2. Steady-state macro aggregates vs. GL Table 4.4¶

The macro aggregates (Y, C, T, YD) arise from the GL Chapter 4 equation system:

$$Y = C + G \qquad T = \theta Y \qquad YD = Y - T + r B_h$$ $$C = \alpha_1 YD + \alpha_2 V \qquad V = V_{-1} + YD - C$$

These are identical to GL regardless of β — only the portfolio split $H_h/B_h$ depends on β. The slight discrepancy in Y arises because the TIR portfolio allocates more wealth to bills (higher yield), which feeds back through $r B_h$ into YD.

In [ ]:

Copied!





model_gl = ModelPC(PCParameters(beta=beta_gl), G_bar=20.0)
ss = model_gl.steady_state(T=300)
V_ss = ss.Hh + ss.Bh

print(f"{'Variable':8s}  {'Simulated':>10s}  {'GL Table 4.4':>12s}  {'Error %':>8s}")
print(f"  {'-'*44}")
for var, gl_val in [("Y",  GL_STEADY_STATE["Y"]),
                    ("C",  GL_STEADY_STATE["C"]),
                    ("T",  GL_STEADY_STATE["T"]),
                    ("YD", GL_STEADY_STATE["YD"])]:
    sim_val = getattr(ss, var)
    err = abs(sim_val - gl_val) / abs(gl_val) * 100
    flag = "✓" if err < 2.0 else "~"
    print(f"  {var:8s}  {sim_val:10.3f}  {gl_val:12.3f}  {err:7.2f}%  {flag}")

print(f"\nPortfolio:")
print(f"  Hh/V = {ss.Hh/V_ss:.4f}  (target {GL_MONEY_SHARE:.4f})")
bs = model_gl.balance_sheet(ss)
print(f"  Balance sheet consistent: {bs.is_consistent()}")
model_gl = ModelPC(PCParameters(beta=beta_gl), G_bar=20.0)
ss = model_gl.steady_state(T=300)
V_ss = ss.Hh + ss.Bh

print(f"{'Variable':8s}  {'Simulated':>10s}  {'GL Table 4.4':>12s}  {'Error %':>8s}")
print(f"  {'-'*44}")
for var, gl_val in [("Y",  GL_STEADY_STATE["Y"]),
                    ("C",  GL_STEADY_STATE["C"]),
                    ("T",  GL_STEADY_STATE["T"]),
                    ("YD", GL_STEADY_STATE["YD"])]:
    sim_val = getattr(ss, var)
    err = abs(sim_val - gl_val) / abs(gl_val) * 100
    flag = "✓" if err < 2.0 else "~"
    print(f"  {var:8s}  {sim_val:10.3f}  {gl_val:12.3f}  {err:7.2f}%  {flag}")

print(f"\nPortfolio:")
print(f"  Hh/V = {ss.Hh/V_ss:.4f}  (target {GL_MONEY_SHARE:.4f})")
bs = model_gl.balance_sheet(ss)
print(f"  Balance sheet consistent: {bs.is_consistent()}")

In [ ]:

Copied!





vars_plot  = ["Y", "C", "T", "YD"]
gl_vals    = [GL_STEADY_STATE[v] for v in vars_plot]
sim_vals   = [getattr(ss, v) for v in vars_plot]

x = np.arange(len(vars_plot))
w = 0.35
fig, ax = plt.subplots(figsize=(7, 4))
ax.bar(x - w/2, gl_vals,  w, label="GL Table 4.4",  color="steelblue",  alpha=0.8)
ax.bar(x + w/2, sim_vals, w, label="EconIAC TIR",   color="darkorange", alpha=0.8)
ax.set_xticks(x); ax.set_xticklabels(vars_plot)
ax.set_ylabel("Value"); ax.set_title("Macro aggregates: GL Table 4.4 vs. EconIAC TIR")
ax.legend()
plt.tight_layout()
plt.show()
vars_plot  = ["Y", "C", "T", "YD"]
gl_vals    = [GL_STEADY_STATE[v] for v in vars_plot]
sim_vals   = [getattr(ss, v) for v in vars_plot]

x = np.arange(len(vars_plot))
w = 0.35
fig, ax = plt.subplots(figsize=(7, 4))
ax.bar(x - w/2, gl_vals,  w, label="GL Table 4.4",  color="steelblue",  alpha=0.8)
ax.bar(x + w/2, sim_vals, w, label="EconIAC TIR",   color="darkorange", alpha=0.8)
ax.set_xticks(x); ax.set_xticklabels(vars_plot)
ax.set_ylabel("Value"); ax.set_title("Macro aggregates: GL Table 4.4 vs. EconIAC TIR")
ax.legend()
plt.tight_layout()
plt.show()

3. Convergence to steady state¶

Model PC has a unique stable steady state for any β ≥ 0. Starting from zero stocks, the economy converges in roughly 60–80 periods. The time path is monotone — no oscillations, no Minsky cycles. (Debt dynamics require the Keen extension.)

In [ ]:

Copied!





traj = model_gl.simulate(T=100)
ts   = range(len(traj))

Y_t  = [s.Y   for s in traj]
Hh_t = [s.Hh  for s in traj]
Bh_t = [s.Bh  for s in traj]
V_t  = [s.Hh + s.Bh for s in traj]

fig, axes = plt.subplots(1, 3, figsize=(13, 4))

axes[0].plot(ts, Y_t, color="steelblue", lw=2)
axes[0].axhline(GL_STEADY_STATE["Y"], ls="--", color="grey", lw=1, label="GL ss")
axes[0].set_xlabel("Period"); axes[0].set_ylabel("Y"); axes[0].set_title("GDP convergence")
axes[0].legend(fontsize=9)

axes[1].plot(ts, Hh_t, color="darkorange", lw=2, label="Hh (money)")
axes[1].plot(ts, Bh_t, color="seagreen",   lw=2, label="Bh (bills)")
axes[1].axhline(GL_STEADY_STATE["Hh"], ls="--", color="darkorange", lw=1, alpha=0.5)
axes[1].axhline(GL_STEADY_STATE["Bh"], ls="--", color="seagreen",   lw=1, alpha=0.5)
axes[1].set_xlabel("Period"); axes[1].set_ylabel("Holdings"); axes[1].set_title("Portfolio convergence")
axes[1].legend(fontsize=9)

# Money share over time
shares_t = [h/(h+b) if (h+b) > 0 else 0.5 for h, b in zip(Hh_t, Bh_t)]
axes[2].plot(ts, shares_t, color="crimson", lw=2)
axes[2].axhline(GL_MONEY_SHARE, ls="--", color="grey", lw=1, label=f"GL share {GL_MONEY_SHARE:.3f}")
axes[2].set_xlabel("Period"); axes[2].set_ylabel("Hh / V"); axes[2].set_title("Money share Hh/V")
axes[2].legend(fontsize=9)

fig.suptitle(f"Model PC convergence  (β = {beta_gl:.1f})", fontsize=13, y=1.02)
plt.tight_layout()
plt.show()
traj = model_gl.simulate(T=100)
ts   = range(len(traj))

Y_t  = [s.Y   for s in traj]
Hh_t = [s.Hh  for s in traj]
Bh_t = [s.Bh  for s in traj]
V_t  = [s.Hh + s.Bh for s in traj]

fig, axes = plt.subplots(1, 3, figsize=(13, 4))

axes[0].plot(ts, Y_t, color="steelblue", lw=2)
axes[0].axhline(GL_STEADY_STATE["Y"], ls="--", color="grey", lw=1, label="GL ss")
axes[0].set_xlabel("Period"); axes[0].set_ylabel("Y"); axes[0].set_title("GDP convergence")
axes[0].legend(fontsize=9)

axes[1].plot(ts, Hh_t, color="darkorange", lw=2, label="Hh (money)")
axes[1].plot(ts, Bh_t, color="seagreen",   lw=2, label="Bh (bills)")
axes[1].axhline(GL_STEADY_STATE["Hh"], ls="--", color="darkorange", lw=1, alpha=0.5)
axes[1].axhline(GL_STEADY_STATE["Bh"], ls="--", color="seagreen",   lw=1, alpha=0.5)
axes[1].set_xlabel("Period"); axes[1].set_ylabel("Holdings"); axes[1].set_title("Portfolio convergence")
axes[1].legend(fontsize=9)

# Money share over time
shares_t = [h/(h+b) if (h+b) > 0 else 0.5 for h, b in zip(Hh_t, Bh_t)]
axes[2].plot(ts, shares_t, color="crimson", lw=2)
axes[2].axhline(GL_MONEY_SHARE, ls="--", color="grey", lw=1, label=f"GL share {GL_MONEY_SHARE:.3f}")
axes[2].set_xlabel("Period"); axes[2].set_ylabel("Hh / V"); axes[2].set_title("Money share Hh/V")
axes[2].legend(fontsize=9)

fig.suptitle(f"Model PC convergence  (β = {beta_gl:.1f})", fontsize=13, y=1.02)
plt.tight_layout()
plt.show()

The function $s(\beta) = \text{gibbs\_weights}([0, r],\, \beta)[0]$ is the calibration landscape: it maps every β to a predicted money share. The curve is:

Monotonically decreasing from 0.5 (at β=0) toward 0 (at β→∞)
Evaluated analytically, not by simulation
One-sided: it cannot reach s > 0.5 at any positive β

This last point is fundamental: it means the TIR model identifies a regime boundary at s = 0.5 that the GL Tobin parameterisation cannot express.

In [ ]:

Copied!





beta_range, shares = portfolio_share_curve(r_bar=0.025)

# Economy observations
economies = {
    "GL (Table 4.4)": (GL_MONEY_SHARE, beta_gl, "steelblue",  "o"),
    "US 2019":        (0.38,           None,    "darkorange", "s"),
    "EA 2019":        (0.45,           None,    "seagreen",   "^"),
    "UK 2019":        (0.52,           None,    "crimson",    "D"),
}

# Calibrate the ones without β
for name, (s_obs, b, *_) in economies.items():
    if b is None:
        b_star, _ = calibrate_beta(s_obs, r_bar=0.025, n_steps=2000, lr=2.0, beta_init=10.0)
        economies[name] = (s_obs, b_star, economies[name][2], economies[name][3])

fig, ax = plt.subplots(figsize=(9, 5))
ax.plot(beta_range, shares, color="navy", lw=2, label="s(β) = Gibbs money share")
ax.axhline(0.5, ls=":", color="grey", lw=1, label="β=0 maximum-entropy (s=0.5)")

for name, (s_obs, b_star, color, marker) in economies.items():
    ax.scatter([b_star], [s_obs], color=color, marker=marker, s=80, zorder=5,
               label=f"{name}  (β*={b_star:.1f}, s={s_obs:.3f})")
    ax.annotate(name, (b_star, s_obs), textcoords="offset points",
                xytext=(6, 4), fontsize=9, color=color)

# Shade the unreachable region s > 0.5
ax.axhspan(0.5, 0.55, alpha=0.08, color="crimson", label="s > 0.5: unreachable at β≥0")

ax.set_xlabel("β (rationality)"); ax.set_ylabel("Money share s = Hh/V")
ax.set_title("Portfolio share curve s(β) at r = 2.5%\nEach economy's β* = implied rationality from national accounts")
ax.set_ylim(0, 0.56); ax.set_xlim(-1, 105)
ax.legend(fontsize=9, loc="upper right")
plt.tight_layout()
plt.show()
beta_range, shares = portfolio_share_curve(r_bar=0.025)

# Economy observations
economies = {
    "GL (Table 4.4)": (GL_MONEY_SHARE, beta_gl, "steelblue",  "o"),
    "US 2019":        (0.38,           None,    "darkorange", "s"),
    "EA 2019":        (0.45,           None,    "seagreen",   "^"),
    "UK 2019":        (0.52,           None,    "crimson",    "D"),
}

# Calibrate the ones without β
for name, (s_obs, b, *_) in economies.items():
    if b is None:
        b_star, _ = calibrate_beta(s_obs, r_bar=0.025, n_steps=2000, lr=2.0, beta_init=10.0)
        economies[name] = (s_obs, b_star, economies[name][2], economies[name][3])

fig, ax = plt.subplots(figsize=(9, 5))
ax.plot(beta_range, shares, color="navy", lw=2, label="s(β) = Gibbs money share")
ax.axhline(0.5, ls=":", color="grey", lw=1, label="β=0 maximum-entropy (s=0.5)")

for name, (s_obs, b_star, color, marker) in economies.items():
    ax.scatter([b_star], [s_obs], color=color, marker=marker, s=80, zorder=5,
               label=f"{name}  (β*={b_star:.1f}, s={s_obs:.3f})")
    ax.annotate(name, (b_star, s_obs), textcoords="offset points",
                xytext=(6, 4), fontsize=9, color=color)

# Shade the unreachable region s > 0.5
ax.axhspan(0.5, 0.55, alpha=0.08, color="crimson", label="s > 0.5: unreachable at β≥0")

ax.set_xlabel("β (rationality)"); ax.set_ylabel("Money share s = Hh/V")
ax.set_title("Portfolio share curve s(β) at r = 2.5%\nEach economy's β* = implied rationality from national accounts")
ax.set_ylim(0, 0.56); ax.set_xlim(-1, 105)
ax.legend(fontsize=9, loc="upper right")
plt.tight_layout()
plt.show()

5. β* calibration across economies¶

Calibrating β to observed Flow of Funds money shares recovers the implied rationality of each household sector — how yield-responsive households are at that economy's interest rate.

Economy	Observed s	β*	Interpretation
GL benchmark	0.270	≈ 40	Moderately rational; high bill allocation
US 2019	0.380	≈ 20	More cash-holding; lower effective β
EA 2019	0.450	≈ 8	Low yield-responsiveness
UK 2019	0.520	≈ 0	Anomaly: s > 0.5 — see below

The novel result is that β* is computable from public data — no survey needed. Flow of Funds tables (Fed Z.1, BoE FA, ECB HFCS) directly identify the TIR rationality parameter for each household sector.

In [ ]:

Copied!





observed = {
    "GL benchmark":         GL_MONEY_SHARE,
    "US households (2019)": 0.38,
    "EA households (2019)": 0.45,
    "UK households (2019)": 0.52,
}

print(f"  {'Economy':28s}  {'Observed':>9s}  {'β*':>8s}  {'Predicted':>10s}  {'Error':>8s}")
print(f"  {'-'*70}")
results = {}
for name, s_obs in observed.items():
    b_star, _ = calibrate_beta(s_obs, r_bar=0.025, n_steps=2000, lr=2.0, beta_init=10.0)
    w_pred    = gibbs_weights(jnp.array([0.0, 0.025]), beta=b_star)
    s_pred    = float(w_pred[0])
    err       = abs(s_pred - s_obs) * 100
    results[name] = (s_obs, b_star, s_pred)
    print(f"  {name:28s}  {s_obs:9.3f}  {b_star:8.2f}  {s_pred:10.4f}  {err:7.4f}%")
observed = {
    "GL benchmark":         GL_MONEY_SHARE,
    "US households (2019)": 0.38,
    "EA households (2019)": 0.45,
    "UK households (2019)": 0.52,
}

print(f"  {'Economy':28s}  {'Observed':>9s}  {'β*':>8s}  {'Predicted':>10s}  {'Error':>8s}")
print(f"  {'-'*70}")
results = {}
for name, s_obs in observed.items():
    b_star, _ = calibrate_beta(s_obs, r_bar=0.025, n_steps=2000, lr=2.0, beta_init=10.0)
    w_pred    = gibbs_weights(jnp.array([0.0, 0.025]), beta=b_star)
    s_pred    = float(w_pred[0])
    err       = abs(s_pred - s_obs) * 100
    results[name] = (s_obs, b_star, s_pred)
    print(f"  {name:28s}  {s_obs:9.3f}  {b_star:8.2f}  {s_pred:10.4f}  {err:7.4f}%")

In [ ]:

Copied!





names   = list(results.keys())
s_obs   = [v[0] for v in results.values()]
b_stars = [v[1] for v in results.values()]
colors  = ["steelblue", "darkorange", "seagreen", "crimson"]

fig, axes = plt.subplots(1, 2, figsize=(12, 4))

axes[0].bar(range(len(names)), s_obs, color=colors, alpha=0.85)
axes[0].axhline(0.5, ls="--", color="grey", lw=1, label="β=0 boundary (s=0.5)")
axes[0].set_xticks(range(len(names))); axes[0].set_xticklabels(names, rotation=15, ha="right")
axes[0].set_ylabel("Observed money share s"); axes[0].set_title("Observed portfolio shares")
axes[0].legend(fontsize=9)

axes[1].bar(range(len(names)), b_stars, color=colors, alpha=0.85)
axes[1].set_xticks(range(len(names))); axes[1].set_xticklabels(names, rotation=15, ha="right")
axes[1].set_ylabel("Implied rationality β*"); axes[1].set_title("Implied rationality β* by economy")

fig.suptitle("TIR β-calibration to Flow of Funds data (r = 2.5%)", fontsize=13)
plt.tight_layout()
plt.show()
names   = list(results.keys())
s_obs   = [v[0] for v in results.values()]
b_stars = [v[1] for v in results.values()]
colors  = ["steelblue", "darkorange", "seagreen", "crimson"]

fig, axes = plt.subplots(1, 2, figsize=(12, 4))

axes[0].bar(range(len(names)), s_obs, color=colors, alpha=0.85)
axes[0].axhline(0.5, ls="--", color="grey", lw=1, label="β=0 boundary (s=0.5)")
axes[0].set_xticks(range(len(names))); axes[0].set_xticklabels(names, rotation=15, ha="right")
axes[0].set_ylabel("Observed money share s"); axes[0].set_title("Observed portfolio shares")
axes[0].legend(fontsize=9)

axes[1].bar(range(len(names)), b_stars, color=colors, alpha=0.85)
axes[1].set_xticks(range(len(names))); axes[1].set_xticklabels(names, rotation=15, ha="right")
axes[1].set_ylabel("Implied rationality β*"); axes[1].set_title("Implied rationality β* by economy")

fig.suptitle("TIR β-calibration to Flow of Funds data (r = 2.5%)", fontsize=13)
plt.tight_layout()
plt.show()

6. The UK anomaly: s > 0.5 and what it means¶

The UK 2019 money share of 0.52 exceeds 0.5 — the maximum-entropy limit of the TIR model at β = 0. This is not a calibration failure; it is a structural finding:

No positive β can produce s > 0.5. The maximum-entropy portfolio is 50/50 at β = 0. Higher β only moves households toward bills, never toward money.

When s > 0.5, households are actively preferring money over the interest-bearing asset, which requires β < 0. This regime has a clear economic interpretation:

Condition	Mechanism
Negative real rates	Bills yield less than inflation; money preferred
Precautionary demand	Liquidity premium exceeds r
Post-QE distortion	Compressed yields remove the incentive to hold bills

The GL Tobin parameterisation cannot express this — lambda0 is bounded by accounting identities, not by a regime boundary. TIR makes the boundary explicit and computable.

This is the novel contribution: the β = 0 line is a phase boundary between yield-seeking (β > 0) and liquidity-seeking (β < 0) household regimes.

In [ ]:

Copied!





# Extend β to negative values to show the liquidity-seeking regime
beta_neg = jnp.linspace(-20.0, 100.0, 600)
U = jnp.array([0.0, 0.025])
shares_ext = jnp.array([float(gibbs_weights(U, beta=float(b))[0]) for b in beta_neg])

fig, ax = plt.subplots(figsize=(10, 5))

# Shade regimes
ax.axvspan(-20, 0,   alpha=0.08, color="crimson",   label="β < 0: liquidity-seeking")
ax.axvspan(0,  100,  alpha=0.05, color="steelblue", label="β > 0: yield-seeking")
ax.axvline(0, color="grey", lw=1.5, ls="--")
ax.axhline(0.5, color="grey", lw=1, ls=":")

ax.plot(beta_neg, shares_ext, color="navy", lw=2.5, label="s(β) extended")

# Mark economies
for name, s_obs, b_star, color, marker in [
    ("GL",  GL_MONEY_SHARE, beta_gl,    "steelblue",  "o"),
    ("US",  0.38,            results["US households (2019)"][1], "darkorange", "s"),
    ("EA",  0.45,            results["EA households (2019)"][1], "seagreen",   "^"),
    ("UK",  0.52,            results["UK households (2019)"][1], "crimson",    "D"),
]:
    ax.scatter([b_star], [s_obs], color=color, marker=marker, s=100, zorder=5)
    ax.annotate(name, (b_star, s_obs), textcoords="offset points",
                xytext=(6, 4), fontsize=10, color=color, fontweight="bold")

ax.set_xlabel("β (rationality, negative = liquidity preference)")
ax.set_ylabel("Money share s = Hh/V")
ax.set_title("Extended portfolio share curve\nβ=0 is a phase boundary between yield-seeking and liquidity-seeking regimes")
ax.set_ylim(0.0, 1.05); ax.set_xlim(-20, 100)
ax.legend(fontsize=9)
plt.tight_layout()
plt.show()

# Show where UK sits relative to the curve
uk_s = 0.52
# Find β that gives closest match (will be near 0 or negative)
idx_closest = int(jnp.argmin(jnp.abs(shares_ext - uk_s)))
print(f"UK s=0.52: closest β on extended curve = {float(beta_neg[idx_closest]):.2f}")
print(f"At β=0: s = {float(gibbs_weights(U, 0.0)[0]):.4f}  (maximum entropy)")
print(f"At β=-5: s = {float(gibbs_weights(U, -5.0)[0]):.4f}  (mild liquidity preference)")
print(f"At β=-15: s = {float(gibbs_weights(U, -15.0)[0]):.4f} (strong liquidity preference)")
# Extend β to negative values to show the liquidity-seeking regime
beta_neg = jnp.linspace(-20.0, 100.0, 600)
U = jnp.array([0.0, 0.025])
shares_ext = jnp.array([float(gibbs_weights(U, beta=float(b))[0]) for b in beta_neg])

fig, ax = plt.subplots(figsize=(10, 5))

# Shade regimes
ax.axvspan(-20, 0,   alpha=0.08, color="crimson",   label="β < 0: liquidity-seeking")
ax.axvspan(0,  100,  alpha=0.05, color="steelblue", label="β > 0: yield-seeking")
ax.axvline(0, color="grey", lw=1.5, ls="--")
ax.axhline(0.5, color="grey", lw=1, ls=":")

ax.plot(beta_neg, shares_ext, color="navy", lw=2.5, label="s(β) extended")

# Mark economies
for name, s_obs, b_star, color, marker in [
    ("GL",  GL_MONEY_SHARE, beta_gl,    "steelblue",  "o"),
    ("US",  0.38,            results["US households (2019)"][1], "darkorange", "s"),
    ("EA",  0.45,            results["EA households (2019)"][1], "seagreen",   "^"),
    ("UK",  0.52,            results["UK households (2019)"][1], "crimson",    "D"),
]:
    ax.scatter([b_star], [s_obs], color=color, marker=marker, s=100, zorder=5)
    ax.annotate(name, (b_star, s_obs), textcoords="offset points",
                xytext=(6, 4), fontsize=10, color=color, fontweight="bold")

ax.set_xlabel("β (rationality, negative = liquidity preference)")
ax.set_ylabel("Money share s = Hh/V")
ax.set_title("Extended portfolio share curve\nβ=0 is a phase boundary between yield-seeking and liquidity-seeking regimes")
ax.set_ylim(0.0, 1.05); ax.set_xlim(-20, 100)
ax.legend(fontsize=9)
plt.tight_layout()
plt.show()

# Show where UK sits relative to the curve
uk_s = 0.52
# Find β that gives closest match (will be near 0 or negative)
idx_closest = int(jnp.argmin(jnp.abs(shares_ext - uk_s)))
print(f"UK s=0.52: closest β on extended curve = {float(beta_neg[idx_closest]):.2f}")
print(f"At β=0: s = {float(gibbs_weights(U, 0.0)[0]):.4f}  (maximum entropy)")
print(f"At β=-5: s = {float(gibbs_weights(U, -5.0)[0]):.4f}  (mild liquidity preference)")
print(f"At β=-15: s = {float(gibbs_weights(U, -15.0)[0]):.4f} (strong liquidity preference)")

7. β-sweep: macro and portfolio response¶

How do the macro aggregates respond as we sweep β? The short answer: Y rises slightly because higher β shifts wealth to bills, increasing interest income $r B_h$, which feeds into YD and hence consumption. But the effect is small — the main action is in the portfolio split.

In [ ]:

Copied!





betas  = [0.0, 1.0, 5.0, 10.0, 20.0, 40.0, beta_gl, 80.0, 120.0]
betas_unique = sorted(set(round(b, 2) for b in betas))

sweep = []
for b in betas_unique:
    m  = ModelPC(PCParameters(beta=b), G_bar=20.0)
    s  = m.steady_state(T=300)
    V  = s.Hh + s.Bh
    ms = s.Hh / V if V > 0 else 0.5
    sweep.append({"beta": b, "Y": s.Y, "Hh": s.Hh, "Bh": s.Bh, "ms": ms})

beta_vals  = [r["beta"] for r in sweep]
Y_vals     = [r["Y"]    for r in sweep]
ms_vals    = [r["ms"]   for r in sweep]
Hh_vals    = [r["Hh"]   for r in sweep]
Bh_vals    = [r["Bh"]   for r in sweep]

fig, axes = plt.subplots(1, 2, figsize=(12, 4))

axes[0].plot(beta_vals, ms_vals, "o-", color="crimson",   lw=2, label="money share Hh/V")
axes[0].axhline(0.5,          ls=":", color="grey",      lw=1)
axes[0].axhline(GL_MONEY_SHARE, ls="--", color="steelblue", lw=1,
                label=f"GL target {GL_MONEY_SHARE:.3f}")
axes[0].set_xlabel("β"); axes[0].set_ylabel("Money share")
axes[0].set_title("Portfolio split vs. β"); axes[0].legend(fontsize=9)

axes[1].plot(beta_vals, Y_vals, "o-", color="steelblue", lw=2, label="Y (GDP)")
axes[1].axhline(GL_STEADY_STATE["Y"], ls="--", color="grey", lw=1, label="GL Y=100")
axes[1].set_xlabel("β"); axes[1].set_ylabel("Y (steady state)")
axes[1].set_title("GDP vs. β\n(small multiplier effect via r·Bh → YD)"); axes[1].legend(fontsize=9)

fig.suptitle("Model PC β-sweep: steady-state sensitivity", fontsize=13)
plt.tight_layout()
plt.show()

# Print table
print(f"  {'β':>8}  {'money share':>12}  {'Y':>8}  {'Hh':>8}  {'Bh':>8}")
print(f"  {'-'*50}")
for r in sweep:
    print(f"  {r['beta']:8.2f}  {r['ms']:12.4f}  {r['Y']:8.3f}  {r['Hh']:8.3f}  {r['Bh']:8.3f}")
betas  = [0.0, 1.0, 5.0, 10.0, 20.0, 40.0, beta_gl, 80.0, 120.0]
betas_unique = sorted(set(round(b, 2) for b in betas))

sweep = []
for b in betas_unique:
    m  = ModelPC(PCParameters(beta=b), G_bar=20.0)
    s  = m.steady_state(T=300)
    V  = s.Hh + s.Bh
    ms = s.Hh / V if V > 0 else 0.5
    sweep.append({"beta": b, "Y": s.Y, "Hh": s.Hh, "Bh": s.Bh, "ms": ms})

beta_vals  = [r["beta"] for r in sweep]
Y_vals     = [r["Y"]    for r in sweep]
ms_vals    = [r["ms"]   for r in sweep]
Hh_vals    = [r["Hh"]   for r in sweep]
Bh_vals    = [r["Bh"]   for r in sweep]

fig, axes = plt.subplots(1, 2, figsize=(12, 4))

axes[0].plot(beta_vals, ms_vals, "o-", color="crimson",   lw=2, label="money share Hh/V")
axes[0].axhline(0.5,          ls=":", color="grey",      lw=1)
axes[0].axhline(GL_MONEY_SHARE, ls="--", color="steelblue", lw=1,
                label=f"GL target {GL_MONEY_SHARE:.3f}")
axes[0].set_xlabel("β"); axes[0].set_ylabel("Money share")
axes[0].set_title("Portfolio split vs. β"); axes[0].legend(fontsize=9)

axes[1].plot(beta_vals, Y_vals, "o-", color="steelblue", lw=2, label="Y (GDP)")
axes[1].axhline(GL_STEADY_STATE["Y"], ls="--", color="grey", lw=1, label="GL Y=100")
axes[1].set_xlabel("β"); axes[1].set_ylabel("Y (steady state)")
axes[1].set_title("GDP vs. β\n(small multiplier effect via r·Bh → YD)"); axes[1].legend(fontsize=9)

fig.suptitle("Model PC β-sweep: steady-state sensitivity", fontsize=13)
plt.tight_layout()
plt.show()

# Print table
print(f"  {'β':>8}  {'money share':>12}  {'Y':>8}  {'Hh':>8}  {'Bh':>8}")
print(f"  {'-'*50}")
for r in sweep:
    print(f"  {r['beta']:8.2f}  {r['ms']:12.4f}  {r['Y']:8.3f}  {r['Hh']:8.3f}  {r['Bh']:8.3f}")

Summary¶

Step	What it showed
β* calibration	GL Table 4.4 money share implies β*≈40 — moderately rational at r=2.5%
Macro replication	Y, C, T, YD agree with GL to ~8% — discrepancy from TIR vs. Tobin feedback
Time series	Monotone convergence in ~60 periods; no Minsky cycles
s(β) curve	Monotone decreasing from 0.5; β=0 is the maximum-entropy limit
Cross-economy β*	US β≈20, EA β≈8, UK β*≈0 (near max-entropy)
UK anomaly	s=0.52 > 0.5 is unreachable at β≥0; identifies a liquidity-preference regime
β-sweep	Small GDP multiplier via r·Bh; main variation is in portfolio split

The key novel result: TIR makes the β=0 phase boundary explicit and computable from public Flow of Funds data. GL's Tobin λ parameters cannot identify this boundary — they are purely descriptive fits with no structural interpretation.

EconIAC is open source: pip install econiac
Theory: Buckley (2026) TIR. doi:10.5281/zenodo.20237288

Tutorial: Godley-Lavoie Model PC — TIR Portfolio Choice¶

1. Calibrate β to GL Table 4.4¶

2. Steady-state macro aggregates vs. GL Table 4.4¶

3. Convergence to steady state¶

4. The portfolio share curve s(β)¶

5. β* calibration across economies¶

6. The UK anomaly: s > 0.5 and what it means¶

7. β-sweep: macro and portfolio response¶

Summary¶