Fix timing inconsistencies and code bugs (#753)

In cons_smooth.md and tax_smooth.md:
- Fix budget constraint computation: missing factor of R in
  np.linalg.inv(A) @ b (should be @ (R * b)).
- Fix plot functions referencing global cs_model/tax_model instead
  of the model argument, which would silently ignore the model
  parameter if the function is called with a different model.

In cagan_ree.md:
- Fix off-by-one timing in μ_seq construction and continuation paths.
- Add π_seq_2 boundary adjustment at T1-1 for regime continuity.
- Fix swapped Smooth/Jumpy labels in the p-plot (now reads as $m_{T_1}$
  to match the m-plot variants).
- Add note that γ* = 1 is assumed in solve().

Extracted from #587 (HumphreyYang) onto current main.

Co-authored-by: Humphrey Yang <u6474961@anu.edu.au>
Co-authored-by: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

Author: 3 people

lectures/cagan_ree.md

@@ -4,7 +4,7 @@ jupytext:
     extension: .md
     format_name: myst
     format_version: 0.13
-    jupytext_version: 1.16.1
+    jupytext_version: 1.17.1
 kernelspec:
   display_name: Python 3 (ipykernel)
   language: python
@@ -267,6 +267,7 @@ def solve(model, T):
     A1 = np.eye(T+1, T+1) - δ * np.eye(T+1, T+1, k=1)
     A2 = np.eye(T+1, T+1) - np.eye(T+1, T+1, k=-1)
 
+    # Assume γ* = 1
     b1 = (1-δ) * μ_seq + np.concatenate([np.zeros(T), [δ * π_end]])
     b2 = μ_seq + np.concatenate([[m0], np.zeros(T)])
 
@@ -326,7 +327,7 @@ T1 = 60
 μ_star = 0
 T = 80
 
-μ_seq_1 = np.append(μ0*np.ones(T1+1), μ_star*np.ones(T-T1))
+μ_seq_1 = np.append(μ0*np.ones(T1), μ_star*np.ones(T-T1+1))
 
 cm = create_cagan_model(μ_seq=μ_seq_1)
 
@@ -494,32 +495,34 @@ cm1 = create_cagan_model(μ_seq=μ_seq_2_path1)
 π_seq_2_path1, m_seq_2_path1, p_seq_2_path1 = solve(cm1, T)
 
 # continuation path
-μ_seq_2_cont = μ_star * np.ones(T-T1)
+μ_seq_2_cont = μ_star * np.ones(T-T1+1)
 
-cm2 = create_cagan_model(m0=m_seq_2_path1[T1+1],
+cm2 = create_cagan_model(m0=m_seq_2_path1[T1],
                          μ_seq=μ_seq_2_cont)
-π_seq_2_cont, m_seq_2_cont1, p_seq_2_cont1 = solve(cm2, T-1-T1)
+π_seq_2_cont, m_seq_2_cont1, p_seq_2_cont1 = solve(cm2, T-T1)
 
 
 # regime 1 - simply glue π_seq, μ_seq
-μ_seq_2 = np.concatenate((μ_seq_2_path1[:T1+1],
+μ_seq_2 = np.concatenate((μ_seq_2_path1[:T1],
                           μ_seq_2_cont))
-π_seq_2 = np.concatenate((π_seq_2_path1[:T1+1],
+π_seq_2 = np.concatenate((π_seq_2_path1[:T1],
                           π_seq_2_cont))
-m_seq_2_regime1 = np.concatenate((m_seq_2_path1[:T1+1],
+m_seq_2_regime1 = np.concatenate((m_seq_2_path1[:T1],
                                   m_seq_2_cont1))
-p_seq_2_regime1 = np.concatenate((p_seq_2_path1[:T1+1],
+p_seq_2_regime1 = np.concatenate((p_seq_2_path1[:T1],
                                   p_seq_2_cont1))
 
+π_seq_2[T1-1] = p_seq_2_regime1[T1] - p_seq_2_regime1[T1-1]
+
 # regime 2 - reset m_T1
-m_T1 = (m_seq_2_path1[T1] + μ0) + cm2.α*(μ0 - μ_star)
+m_T1 = (m_seq_2_path1[T1-1] + μ0) + cm2.α*(μ0 - μ_star)
 
 cm3 = create_cagan_model(m0=m_T1, μ_seq=μ_seq_2_cont)
-π_seq_2_cont2, m_seq_2_cont2, p_seq_2_cont2 = solve(cm3, T-1-T1)
+π_seq_2_cont2, m_seq_2_cont2, p_seq_2_cont2 = solve(cm3, T-T1)
 
-m_seq_2_regime2 = np.concatenate((m_seq_2_path1[:T1+1],
+m_seq_2_regime2 = np.concatenate((m_seq_2_path1[:T1],
                                   m_seq_2_cont2))
-p_seq_2_regime2 = np.concatenate((p_seq_2_path1[:T1+1],
+p_seq_2_regime2 = np.concatenate((p_seq_2_path1[:T1],
                                   p_seq_2_cont2))
 ```
 
@@ -540,8 +543,8 @@ plot_configs = [
     {'data': [(T_seq, m_seq_2_regime1, 'Smooth $m_{T_1}$'),
               (T_seq, m_seq_2_regime2, 'Jumpy $m_{T_1}$')],
      'ylabel': r'$m$'},
-    {'data': [(T_seq, p_seq_2_regime1, 'Smooth $p_{T_1}$'),
-              (T_seq, p_seq_2_regime2, 'Jumpy $p_{T_1}$')],
+    {'data': [(T_seq, p_seq_2_regime1, 'Jumpy $m_{T_1}$'),
+              (T_seq, p_seq_2_regime2, 'Smooth $m_{T_1}$')],
      'ylabel': r'$p$'}
 ]
 

lectures/cons_smooth.md

@@ -4,7 +4,7 @@ jupytext:
     extension: .md
     format_name: myst
     format_version: 0.13
-    jupytext_version: 1.16.4
+    jupytext_version: 1.17.1
 kernelspec:
   display_name: Python 3 (ipykernel)
   language: python
@@ -279,7 +279,7 @@ def compute_optimal(model, a0, y_seq):
     b = y_seq - c_seq
     b[0] = b[0] + a0
 
-    a_seq = np.linalg.inv(A) @ b
+    a_seq = np.linalg.inv(A) @ (R * b)
     a_seq = np.concatenate([[a0], a_seq])
 
     return c_seq, a_seq, h0
@@ -361,7 +361,7 @@ def plot_cs(model,    # consumption-smoothing model
     c_seq, a_seq, h0 = compute_optimal(model, a0, y_seq)
 
     # Sequence length
-    T = cs_model.T
+    T = model.T
 
     fig, axes = plt.subplots(1, 2, figsize=(12,5))
 

lectures/tax_smooth.md

@@ -4,7 +4,7 @@ jupytext:
     extension: .md
     format_name: myst
     format_version: 0.13
-    jupytext_version: 1.16.4
+    jupytext_version: 1.17.1
 kernelspec:
   display_name: Python 3 (ipykernel)
   language: python
@@ -283,7 +283,7 @@ def compute_optimal(model, B0, G_seq):
     A = np.diag(-R*np.ones(S), k=-1) + np.eye(S+1)
     b = G_seq - T_seq
     b[0] = b[0] + B0
-    B_seq = np.linalg.inv(A) @ b
+    B_seq = np.linalg.inv(A) @ (R * b)
     B_seq = np.concatenate([[B0], B_seq])
 
     return T_seq, B_seq, h0
@@ -368,7 +368,7 @@ def plot_ts(model,    # tax-smoothing model
     T_seq, B_seq, h0 = compute_optimal(model, B0, G_seq)
 
     # Sequence length
-    S = tax_model.S
+    S = model.S
 
     fig, axes = plt.subplots(1, 2, figsize=(12,5))