WIP

Author: albop

ICAI/partial.dyno

@@ -17,6 +17,9 @@ r <- r_w
 η <- 2.5
 η_b <- 2.5
 ω <- 0.5
+
+T <- 20
+
 hand_to_mouth <- 0.0
 
 y[~] <- y_bar
@@ -31,47 +34,59 @@ b_τ[~] <- 0
 p[~] <- β/(1-β)*r
 q[~] <- p[~]*a[~]*y[~]/r
 n[~] <- 1 # total number of shares
-x[~] <- 1/χ # number of shares per top earner
 
-W_τ[~] <- p[~]*(b_τ[~]) + q[~]*x[~]
+W_τ[~] <- p[~]*(b_τ[~]) + q[~]
 W_b[~] <- 1.0 + p[~]*(b_b[~])*0  # adding/removing the last term affects the result !
-c_b[~] <- a[~]*(1-y[~])*(1-ζ[~])/(1-χ) 
-c_τ[~] <- a[~]*(1-y[~])*ζ[~]*x[~] + b_τ[~]*r + (a[~]*y[~])*x[~]
-y_b[~] <- a[~]*(1-y[~])*(1-ζ[~])/(1-χ)    # labor income for bottom-earners
-y_τ[~] <- a[~]*(1-y[~])*ζ[~]*x[~] + (a[~]*y[~])*x[~]     # top earners income (labor plus dividends)
+c_b[~] <- a[~]*(1-y[~])*(1-ζ[~])
+c_τ[~] <- a[~]*(1-y[~])*ζ[~] + b_τ[~]*r + (a[~]*y[~])
+y_b[~] <- a[~]*(1-y[~])*(1-ζ[~])    # labor income for bottom-earners
+y_τ[~] <- a[~]*(1-y[~])*ζ[~] + (a[~]*y[~])     # top earners income (labor plus dividends)
 c_b[~] <- b_b[~]*r + y_b[~]
 c_τ[~] <- b_τ[~]*r + y_τ[~]
 
-
 φ_τ <- (1-β_τ/p[~])*c_τ[~]^(-σ)/(W_τ[~])^(-η)
 φ_b <- (1-β_b/p[~])*c_b[~]^(-σ)/(W_b[~])^(-η_b)
 
+ΔNFA[~] <- 0.0
 
 
 # exogenous processes
 y[t] = y[~] + ρ_y*(y[t-1]-y[~])
 ζ[t] = ζ[~] + ρ_ζ*(ζ[t-1]-ζ[~])
-x[t] = x[~]
+
 a[t] = 1 + ρ_a*(a[t-1]-1)
 
-d[t]   = d[t-1]   + e_d[t]    # income shock for top earners
-d_b[t] = d_b[t-1] + e_d_b[t]  # income shock for bottom earners
+d[t]   = d[t-1]   + e_d[t]    # income shock for top earners (aggregate)
+d_b[t] = d_b[t-1] + e_d_b[t]  # income shock for bottom earners (aggregate)
 
-y_b[t] = a[t]*(1-y[t])*(1-ζ[t])/(1-χ)    # labor income for bottom-earners
-y_τ[t] = a[t]*(1-y[t])*ζ[t]*x[t] + (a[t]*y[t])*x[t]     # top earners income (labor plus dividends)
+y_b[t] = a[t]*(1-y[t])*(1-ζ[t])                      # labor income for bottom-earners (aggregate)
+y_τ[t] = a[t]*(1-y[t])*ζ[t] + (a[t]*y[t])            # top earners income (per capita)
 
-c_b[t] = d_b[t]/(1-χ) +  b_b[t-1]*r - (b_b[t]-b_b[t-1])*p[t] + y_b[t]
-c_τ[t] = d[t]         +  b_τ[t-1]*r - (b_τ[t]-b_τ[t-1])*p[t] + y_τ[t]
+c_b[t] = d_b[t] +  b_b[t-1]*r - (b_b[t]-b_b[t-1])*p[t] + y_b[t]  # consumption bottom earners (aggregate)
+c_τ[t] =   d[t] +  b_τ[t-1]*r - (b_τ[t]-b_τ[t-1])*p[t] + y_τ[t]  # consumption top earners    (aggregate)
 
-W_τ[t] = q[t]*x[t] + p[t]*(b_τ[t])
-W_b[t] = 1         + p[t]*(b_b[t])
+W_τ[t] = q[t] + p[t]*(b_τ[t])     # wealth top earners (aggregate)
+W_b[t] = 1         + p[t]*(b_b[t])     # wealth bottom earners (aggregate)
 
 p[t] = β_b*(c_b[t+1]/c_b[t])^(-σ)*(r+p[t+1]) + φ_b*(W_b[t])^(-η_b)/(c_b[t])^(-σ)*p[t]
 p[t] = β_τ*(c_τ[t+1]/c_τ[t])^(-σ)*(r+p[t+1]) + φ_τ*(W_τ[t])^(-η  )/(c_τ[t])^(-σ)*p[t]
 
+q[t] = p[~]*a[~]*y[~]/r
+
+ΔNFA[t] = b_τ[t] + b_b[t] - (b_τ[~] + b_b[~])   # change in net foreign asset position
+
+# small open economy assumption:
 r_w + p[t] = 1
 
-q[t] = p[~]*a[~]*y[~]/r
+# closed economy assumption:
+# b_τ[t] + b_b[t] = b_τ[~] + b_b[~] # net-zero change in bond supply
+
+
+### Exogenous shocks
+# e_d[t]   <- N(0, 0.01)
+# e_d_b[t] <- N(0, 0.01)
+
 
-# e_d[t]   <- N(0, 0.1)
-e_d_b[t] <- N(0, 0.1)
+e_d[1] <- 0.01
+forall t, 2<=t<T: e_d[t] <- 0.02
+e_d_b[1] <- 0.01

ICAI/solve_icai.py

@@ -1,10 +1,64 @@
 import dyno
 from matplotlib import pyplot as plt
 from dyno.symbolic_model import SymbolicModel
+import pandas as pd
 
 model = SymbolicModel("partial.dyno")
 
+v0 = model.deterministic_guess()
+
+pd.DataFrame(v0, columns=model.symbols['variables'])
+
+
+
+
+model.jacobians[0]
+
+
 dr = model.solve()
+sim = dr.irfs(type='level')['e_d']
+# plt.plot(sim.index, sim['W_τ'], label="baseline")
+
+plt.plot(sim.index, sim['b_b'], label="baseline",marker='.',linestyle='-')
+plt.plot(sim.index, sim['b_τ'], label="baseline",marker='.',linestyle='-')
+plt.xlim(0,10)
+
+plt.plot(sim.index, sim['W_b'], label="baseline",marker='.',linestyle='-')
+plt.plot(sim.index, sim['W_τ'], label="baseline",marker='.',linestyle='-')
+plt.xlim(0,10)
+
+
+m2 = model.recalibrate(η=1.51)
+dr2 = m2.solve()
+sim2 = dr2.irfs()['e_d']
+
+# plt.plot(sim.index, (sim['W_τ']-sim['W_τ'][0])/sim['d'], label="baseline")
+# plt.plot(sim2.index, (sim2['W_τ']-sim2['W_τ'][0])/sim2['d'], label="η=3.0")
+
+plt.plot(sim2.index, sim2['W_τ'])
+
+
+m2 = model.recalibrate(η=1.51)
+dr2 = m2.solve()
+sim2 = dr2.irfs()['e_d']
+
+# plt.plot(sim.index, (sim['W_τ']-sim['W_τ'][0])/sim['d'], label="baseline")
+# plt.plot(sim2.index, (sim2['W_τ']-sim2['W_τ'][0])/sim2['d'], label="η=3.0")
+
+plt.plot(sim.index, (sim['W_b']-sim['W_b'][0])/sim['d'], label="baseline")
+plt.plot(sim2.index, (sim2['W_b']-sim2['W_b'][0])/sim2['d'], label="η=3.0")
+
+
+
+
+
+print(m2.data.context['constants']['D'])
+print(m2.data.context['processes'][('e_d',)].Σ)
+
+
+plt.plot(sim['e_d'].index, sim['e_d']['d'], label="baseline")
+plt.plot(sim2['e_d'].index, sim2['e_d']['d'], label="D=0.05")
+plt.legend()
 
 
 ###### 

dev_dynare.py

@@ -4,16 +4,52 @@
 dir_path = os.path.dirname(os.path.realpath(__file__))
 
 
-f_r = "examples/modfiles/RBC.mod"
 f_o = "examples/RBC.dyno"
+f_r = "examples/modfiles/RBC.mod"
+
+
+from dyno.symbolic_model import SymbolicModel
+
+model1 = SymbolicModel(f_o)
+model2 = SymbolicModel(f_r)
+
+
+import copy
+
+t1 = time.time()
+m1 = model1.recalibrate(beta=0.4)
+t2 = time.time()
+print("Elapsed: ", t2 - t1)
+
+
+
+import copy
+t1 = time.time()
+mm1 = copy.deepcopy(model1)
+mm1.solve()
+t2 = time.time()
+print("Elapsed: ", t2 - t1)
+
+
+import copy
+t1 = time.time()
+mm2 = copy.deepcopy(model2)
+dr = mm2.solve()
+t2 = time.time()
+print("Elapsed: ", t2 - t1)
+
+
 
 from dyno.symbolic_model import SymbolicModel
 
 t1 = time.time()
 model = SymbolicModel(f_o)
 t2 = time.time()
 dr1 = model.solve()
+t3 = time.time()
 print("Model parsing time: ", t2 - t1)
+print("Model solving time: ", t3 - t2)
+
 
 from dyno.symbolic_model import SymbolicModel
 
@@ -23,7 +59,12 @@
 dr2 = model.solve()
 t3 = time.time()
 print("Model parsing time: ", t2 - t1)
-# print("Model solving time: ", t3 - t2)
+print("Model solving time: ", t3 - t2)
+
+
+
+
+
 
 from dyno.modfile import DynareModel as DynareModelPP
 
@@ -33,7 +74,7 @@
 dr3 = model.solve()
 t3 = time.time()
 print("Model parsing time: ", t2 - t1)
-# print("Model solving time: ", t3 - t2)
+print("Model solving time: ", t3 - t2)
 
 
 from IPython.display import Math, display