Step 1 - SIM

1 Shared Setup

Load minimal packages and global options only.

Check that required packages are installed. Stop early with a clear error if anything is missing. This keeps the tutorial run deterministic and minimal.

required_packages <- c("jsonlite", "ggplot2")
missing_packages <- required_packages[!vapply(required_packages, requireNamespace, logical(1), quietly = TRUE)]
if (length(missing_packages) > 0) {
  stop("Missing package(s): ", paste(missing_packages, collapse = ", "), call. = FALSE)
}

1.1 Function and Variable Location Map

This step defines everything locally in this file: steps/step01_sim.qmd.

• simulate_sim_lag() is defined in chunk step1_functions_sim_core.
• simulate_sim_current_fp() is defined in chunk step1_functions_sim_current.
• Baseline SIM parameters (G, theta, alpha1, alpha2, Y0, H0) are set in chunk step1_run_lag_below.
• Comparison run objects (lag, cur) are set in chunk step1_run_compare_lag_current.
• Core play variable (Y0_play) is set in chunk step1_play_core_initial.
• Optional play variable (alpha2_play) is set in chunk step1_play_optional_alpha2.

2 Step 1: SIM Closed Model (Chapter 3 of Godley & Lavoie 2012 book)

2.1 Objective

Build a simple (SIM) SFC model with Government, Production and Household sectors. We see that GDP and Wealth reach a steady state after some time! Depending on whether consumption behaviour depends on current year vs previous year income, we need to solve differently

2.2 Why this step matters

This is the accounting foundation for all later steps. If we clearly see how deficits, wealth accumulation and output adjustment work in SIM, then adding IOT, RoW and emissions later is much easier to understand.

2.3 What changed from previous step

This is the first economic model step, so there is no previous dynamics to compare with. In this step we introduce the baseline accounting logic and two timing variants (SIM-lag and SIM-current).

2.4 Equations

• Output identity Y_t = C_t + G_t
• Wages W_t = Y_t , combining wages and profits as wages to households
• Taxes T_t = theta * Y_t
• Household Disposable Income YD_t = Y_t - T_t
• Consumption C_t(lag) = alpha1 * YD_{t-1} + alpha2 * H_{t-1} , (assume based on last year’s income, SIM-lag)
• Consumption C_t(current) = alpha1 * YD_t + alpha2 * H_{t-1} , (assume based on current year’s income, SIM-current)
• Household Wealth / change in wealth H_t = H_{t-1} + YD_t - C_t
• Government deficit/debt: DEF_t = G_t - T_t, B_t = B_{t-1} + DEF_t
• Fixed point of GDP: Y* = G / theta
• Fixed point consumption: C* = (1-theta) * Y*

Transaction Flow Matrix (TFM, one period, sign by sector column):

Flow	Households	Government	Production	Sum
Wages / income Y	+Y	0	-Y	0
Taxes T	-T	+T	0	0
Consumption C	-C	0	+C	0
Government demand G	0	-G	+G	0
Change in wealth/debt Delta_H (= Delta_B)	-Delta_H	+Delta_B	0	0
Column sum	0	0	0	0

From column sums: - Households: Y - T - C - Delta_H = 0, so Delta_H = Y - T - C. - Government: T - G + Delta_B = 0, so Delta_B = G - T. - Production: -Y + C + G = 0, so Y = C + G. - Stock-flow consistency: Delta_H = Delta_B.

Column-sum check (explicit): household column sum = 0, government column sum = 0, production column sum = 0.

Balance Sheet Matrix (BSM, simplified SIM view):

	Households	Government	Production	Sum
Wealth/debt stock H (=B)	+H	-B	0	0
Net worth	V_h	V_g	V_p	0

Usual BSM convention: each sector column sums to its net worth (not zero). If a balancing net-worth row is written with opposite signs, then each column sums to zero.

Net-worth identities: - V_h = H - V_g = -B - V_p = 0 - with H = B, V_h + V_g + V_p = 0

Steady-state derivation (SIM): - At steady state, stocks are constant: H_t = H_{t-1} = H* so 0 = YD* - C*, hence C* = YD* = (1-theta)Y*. - GDP: Y* = C* + G = (1-theta)Y* + G*. - Rearranging: theta * Y* = G, so Y* = G / theta. - From behavior at steady state: C* = alpha1*(1-theta)Y* + alpha2*H* and C* = (1-theta)Y*. - Therefore: H* = ((1-alpha1)*(1-theta)/alpha2) * Y* = ((1-alpha1)*(1-theta)/alpha2) * (G/theta). - The fixed-point GDP Y* is independent of alpha1 and alpha2; those affect only transition dynamics.

SIM-lag is easy to simulate, since all updates depend on previous year’s values, so no circular depndencies

Define the SIM-lag dynamics as a simple time loop. Consumption C_t(lag) = alpha1 * YD_{t-1} + alpha2 * H_{t-1} depends on lagged disposable income and past wealth. The function simulates the dynamics for time steps T (years), and returns a time-series for plotting and comparisons.

simulate_sim_lag <- function(T = 20, G = 100, theta = 0.2, alpha1 = 0.6, alpha2 = 0.2, Y0 = 80, H0 = 300) {
  out <- data.frame(t = seq_len(T), Y = NA_real_, T_tax = NA_real_, YD = NA_real_, C = NA_real_, H = NA_real_)
  YD_prev <- (1 - theta) * Y0
  H_prev <- H0
  for (tt in seq_len(T)) {
    C_t <- alpha1 * YD_prev + alpha2 * H_prev
    Y_t <- C_t + G
    T_t <- theta * Y_t
    YD_t <- Y_t - T_t
    H_t <- H_prev + YD_t - C_t
    out[tt, c("Y", "T_tax", "YD", "C", "H")] <- c(Y_t, T_t, YD_t, C_t, H_t)
    YD_prev <- YD_t
    H_prev <- H_t
  }
  out
}

Run SIM-lag from initial values that are below the steady state.

Set baseline SIM parameters and compute the fixed point Y. Choose an initial GDP below Y and derive a consistent H0. Run the lag model to generate the convergence path.

T_h <- 40 # time horizon
G <- 100  # government spending per year
theta <- 0.2 # tax rate
alpha1 <- 0.6 # fraction consumption out of income
alpha2 <- 0.2 # fraction consumption out of wealth
# fixed points for GDP and wealth
Y_star <- G / theta
C_star <- (1 - theta) * Y_star
H_star <- ((1 - alpha1) * (1 - theta) * Y_star) / alpha2

# setting Y0 below the fixed point, try above the fixed point as well!
Y0 <- 0.7 * Y_star
#Y0 <- 1.3 * Y_star

# set H0 consistent with Y0, try something else, e.g. 0.0 as well!
H0 <- (Y0 * (1 - alpha1 * (1 - theta)) - G) / alpha2
#H0 <- 0.0

# also define above-fixed-point initial values for comparison plots
lag <- simulate_sim_lag(T = T_h, G = G, theta = theta, alpha1 = alpha1, alpha2 = alpha2, Y0 = Y0, H0 = H0)
lag$DEF <- G - lag$T_tax
lag$H_from_DEF <- H0 + cumsum(lag$DEF)

Plot Step 1A accounting paths in one figure.

This plot overlays GDP, consumption, deficit and wealth identities to show the core stock-flow link clearly: government deficit flow accumulates into household wealth.

y_all <- c(lag$Y, lag$C, lag$H, lag$H_from_DEF, lag$DEF, Y_star, C_star, H_star)
plot(lag$t, lag$Y, type = "l", lwd = 2, col = "steelblue",
     ylim = range(y_all, na.rm = TRUE),
     xlab = "t", ylab = "level", main = "Step 1A: SIM-lag accounting dynamics")
lines(lag$t, lag$C, lwd = 2, col = "darkgreen")
lines(lag$t, lag$DEF, lwd = 2, col = "firebrick")
lines(lag$t, lag$H_from_DEF, lwd = 4, col = "goldenrod")
lines(lag$t, lag$H, lwd = 2, col = "black")
abline(h = Y_star, lty = 2, col = "steelblue")
abline(h = C_star, lty = 2, col = "darkgreen")
abline(h = H_star, lty = 2, col = "black")
legend("right",
       legend = c("GDP Y", "Consumption C", "Deficit (G - Tax)", "H0 + cumulative deficit", "Wealth H", "Y*", "C*", "H*"),
       col = c("steelblue", "darkgreen", "firebrick", "goldenrod", "black", "steelblue", "darkgreen", "black"),
       lty = c(1, 1, 1, 1, 1, 2, 2, 2), lwd = c(2, 2, 2, 4, 2, 1, 1, 1), bty = "n")

Core play: See if same fixed point is reached with a different initial condition, e.g. above the fixed point

How about if we choose H0 that is not “consistent” with the Y0, e.g. 0.0?

SIM-current with Consumption C_t(current) = alpha1 * YD_t + alpha2 * H_{t-1} is not so easy to simulate, since updates are interdependent. C_t depends on current disposable income YD_t and past wealth. Current disposable income YD_t = Y_t - T_t = C_t + G - T_t. Thus, there is a circularity/interdependence. So we need to solve for these variables consistently with each other.

For a simple linear system, we can invert directly the system of equations directly. But in general, we perform fixed-point iteration within each period / time step, here a year. I n this specific linear SIM-current case, the system reduces to Y_t = G + alpha1 * (1 - theta) * Y_t + alpha2 * H_{t-1}. So we can solve directly as Y_t = (G + alpha2 * H_{t-1}) / (1 - alpha1 * (1 - theta)). We still use fixed-point iteration on the full block (C -> Y -> T -> YD -> C) because that approach extends directly to nonlinear or richer model variants.

simulate_sim_current_fp <- function(T = 20, G = 100, theta = 0.2, alpha1 = 0.6, alpha2 = 0.2,
                                    Y0 = 80, H0 = 300, fp_tol = 1e-10, fp_max = 200) {
  out <- data.frame(t = seq_len(T), Y = NA_real_, T_tax = NA_real_, YD = NA_real_, C = NA_real_, H = NA_real_)
  H_prev <- H0
  C_prev <- Y0 - G
  for (tt in seq_len(T)) {
    C_guess <- C_prev
    n_iter <- 0L
    repeat {
      n_iter <- n_iter + 1L
      Y_t <- C_guess + G
      T_t <- theta * Y_t
      YD_t <- Y_t - T_t
      C_new <- alpha1 * YD_t + alpha2 * H_prev
      if (abs(C_new - C_guess) < fp_tol || n_iter >= fp_max) break
      C_guess <- C_new
    }
    C_t <- C_new
    Y_t <- C_t + G
    T_t <- theta * Y_t
    YD_t <- Y_t - T_t
    H_t <- H_prev + YD_t - C_t
    out[tt, c("Y", "T_tax", "YD", "C", "H")] <- c(Y_t, T_t, YD_t, C_t, H_t)
    H_prev <- H_t
    C_prev <- C_t
  }
  out
}

Compare lag and current dynamics.

Run SIM-current with the same initial condition used above. Plot lag and current dynamics together so participants can compare timing assumptions under one shared starting point.

lag <- simulate_sim_lag(T = T_h, G = G, theta = theta, alpha1 = alpha1, alpha2 = alpha2, Y0 = Y0, H0 = H0)
cur <- simulate_sim_current_fp(T = T_h, G = G, theta = theta, alpha1 = alpha1, alpha2 = alpha2, Y0 = Y0, H0 = H0)

yr <- range(c(lag$Y, cur$Y, Y_star), na.rm = TRUE)
pad <- 0.05 * diff(yr); if (!is.finite(pad) || pad <= 0) pad <- 1
plot(lag$t, lag$Y, type = "l", lwd = 2, col = "steelblue",
     ylim = c(yr[1] - pad, yr[2] + pad),
     xlab = "t", ylab = "GDP", main = "Step 1C: SIM-lag vs SIM-current")
lines(cur$t, cur$Y, lwd = 2, col = "darkgreen")
abline(h = Y_star, lty = 2)
legend("right", legend = c("Lag", "Current", "Y*"),
       col = c("steelblue", "darkgreen", "black"),
       lty = c(1, 1, 2), bty = "n")

Change initial GDP and compare convergence across lag and current versions.

Set Y0_play above fixed point or exactly at fixed point, and compare lag/current convergence under the same parameters.

Y0_play <- 1.3 * Y_star  # @exercise[id=step1_y0;kind=core;question_expr=1.3 * Y_star;prompt="Change initial GDP and compare convergence (below/above/equal to fixed point)";hint="Try 0.7*Y_star, 1.0*Y_star, 1.3*Y_star"]
H0_play <- (Y0_play * (1 - alpha1 * (1 - theta)) - G) / alpha2

lag_play <- simulate_sim_lag(T = T_h, G = G, theta = theta, alpha1 = alpha1, alpha2 = alpha2, Y0 = Y0_play, H0 = H0_play)
cur_play <- simulate_sim_current_fp(T = T_h, G = G, theta = theta, alpha1 = alpha1, alpha2 = alpha2, Y0 = Y0_play, H0 = H0_play)

plot(lag_play$t, lag_play$Y, type = "l", lwd = 2, col = "steelblue",
     xlab = "t", ylab = "GDP", main = "Step 1 Core Play: change initial GDP")
lines(cur_play$t, cur_play$Y, lwd = 2, col = "darkgreen")
abline(h = Y_star, lty = 2)
legend("right", legend = c("Lag", "Current", "Y*"),
       col = c("steelblue", "darkgreen", "black"),
       lty = c(1, 1, 2), bty = "n")

Optional play: adjust wealth propensity in both variants.

Vary the propensity to consume out of wealth (alpha2). Rerun both lag and current versions with the new alpha2. Plot the impact on speed of convergence.

Try also alpha2 = 0.0!

alpha2_play <- 0.25  # @exercise[id=step1_alpha2;kind=optional;question_expr=0.1;prompt="Change alpha2 and compare speed of convergence";hint="Lower alpha2 usually slows wealth-driven demand"]

lag_opt <- simulate_sim_lag(T = T_h, G = G, theta = theta, alpha1 = alpha1, alpha2 = alpha2_play, Y0 = Y0, H0 = H0)
cur_opt <- simulate_sim_current_fp(T = T_h, G = G, theta = theta, alpha1 = alpha1, alpha2 = alpha2_play, Y0 = Y0, H0 = H0)

plot(lag_opt$t, lag_opt$Y, type = "l", lwd = 2, col = "steelblue",
     xlab = "t", ylab = "GDP", main = "Step 1 Optional: effect of alpha2")
lines(cur_opt$t, cur_opt$Y, lwd = 2, col = "darkgreen")
legend("right", legend = c("Lag", "Current"), col = c("steelblue", "darkgreen"), lty = 1, bty = "n")

2.5 Insight (Step 1)

• G and theta determine the fixed point location.
• alpha1 and alpha2 mainly shape transition speed and path.
• The deficit flow and wealth accumulation track each other by accounting identity.

2.6 What can go wrong / interpretation caveat

• Inconsistent initial values (Y0, H0) can create strong transients.
• This is a stripped-down pedagogical model: no prices, banking, inventories, or firm financial structure.

2.7 Interpretation

Both dynamic versions converge to the same fixed point, but with different transient paths.