Skip to content

Commit

Permalink
updates
Browse files Browse the repository at this point in the history
  • Loading branch information
@annabellasd
annabellasd committed Oct 10, 2023
1 parent 088e6a2 commit befd388
Show file tree
Hide file tree
Showing 11 changed files with 343 additions and 349 deletions.
50 changes: 25 additions & 25 deletions lectures/multi_agent_models/aiyagari.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -210,23 +210,23 @@ using LaTeXStrings, Parameters, Plots, QuantEcon
```

```{code-cell} julia
Household = @with_kw (r = 0.01,
w = 1.0,
σ = 1.0,
β = 0.96,
z_chain = MarkovChain([0.9 0.1; 0.1 0.9], [0.1; 1.0]),
a_min = 1e-10,
a_max = 18.0,
a_size = 200,
a_vals = range(a_min, a_max, length = a_size),
z_size = length(z_chain.state_values),
n = a_size * z_size,
s_vals = gridmake(a_vals, z_chain.state_values),
s_i_vals = gridmake(1:a_size, 1:z_size),
u = σ == 1 ? x -> log(x) : x -> (x^(1 - σ) - 1) / (1 - σ),
R = setup_R!(fill(-Inf, n, a_size), a_vals, s_vals, r, w, u),
# -Inf is the utility of dying (0 consumption)
Q = setup_Q!(zeros(n, a_size, n), s_i_vals, z_chain))
Household(;r = 0.01,
w = 1.0,
sigma = 1.0,
beta = 0.96,
z_chain = MarkovChain([0.9 0.1; 0.1 0.9], [0.1; 1.0]),
a_min = 1e-10,
a_max = 18.0,
a_size = 200,
a_vals = range(a_min, a_max, length = a_size),
z_size = length(z_chain.state_values),
n = a_size * z_size,
s_vals = gridmake(a_vals, z_chain.state_values),
s_i_vals = gridmake(1:a_size, 1:z_size),
u = sigma == 1 ? x -> log(x) : x -> (x^(1 - sigma) - 1) / (1 - sigma),
R = setup_R!(fill(-Inf, n, a_size), a_vals, s_vals, r, w, u),
# -Inf is the utility of dying (0 consumption)
Q = setup_Q!(zeros(n, a_size, n), s_i_vals, z_chain)) = (; r, w, sigma, beta, z_chain, a_min, a_max, a_size, a_vals, z_size, n, s_vals, s_i_vals, u, R, Q)
function setup_Q!(Q, s_i_vals, z_chain)
for next_s_i in 1:size(Q, 3)
Expand Down Expand Up @@ -274,7 +274,7 @@ Random.seed!(42);
am = Household(a_max = 20.0, r = 0.03, w = 0.956)
# Use the instance to build a discrete dynamic program
am_ddp = DiscreteDP(am.R, am.Q, am.β)
am_ddp = DiscreteDP(am.R, am.Q, am.beta)
# Solve using policy function iteration
results = solve(am_ddp, PFI)
Expand Down Expand Up @@ -326,16 +326,16 @@ Random.seed!(42);
# Firms' parameters
const A = 1
const N = 1
const α = 0.33
const β = 0.96
const δ = 0.05
const alpha = 0.33
const beta = 0.96
const delta = 0.05
function r_to_w(r)
return A * (1 - α) * (A * α / (r + δ)) ^ (α / (1 - α))
return A * (1 - alpha) * (A * alpha / (r + delta)) ^ (alpha / (1 - alpha))
end
function rd(K)
return A * α * (N / K) ^ (1 - α) - δ
return A * alpha * (N / K) ^ (1 - alpha) - delta
end
function prices_to_capital_stock(am, r)
Expand All @@ -345,7 +345,7 @@ function prices_to_capital_stock(am, r)
(;a_vals, s_vals, u) = am
setup_R!(am.R, a_vals, s_vals, r, w, u)
aiyagari_ddp = DiscreteDP(am.R, am.Q, am.β)
aiyagari_ddp = DiscreteDP(am.R, am.Q, am.beta)
# Compute the optimal policy
results = solve(aiyagari_ddp, PFI)
Expand All @@ -358,7 +358,7 @@ function prices_to_capital_stock(am, r)
end
# Create an instance of Household
am = Household(β = β, a_max = 20.0)
am = Household(beta = beta, a_max = 20.0)
# Create a grid of r values at which to compute demand and supply of capital
r_vals = range(0.005, 0.04, length = 20)
Expand Down
69 changes: 34 additions & 35 deletions lectures/multi_agent_models/arellano.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -312,56 +312,55 @@ using Test
```

```{code-cell} julia
function ArellanoEconomy(;β = .953,
γ = 2.,
function ArellanoEconomy(;beta = .953,
gamma = 2.,
r = 0.017,
ρ = 0.945,
η = 0.025,
θ = 0.282,
rho = 0.945,
eta = 0.025,
theta = 0.282,
ny = 21,
nB = 251)
# create grids
Bgrid = collect(range(-.4, .4, length = nB))
mc = tauchen(ny, ρ, η)
Π = mc.p
mc = tauchen(ny, rho, eta)
Pi = mc.p
ygrid = exp.(mc.state_values)
ydefgrid = min.(.969 * mean(ygrid), ygrid)
# define value functions
# notice ordered different than Python to take
# advantage of column major layout of Julia)
# advantage of column major layout of Julia
vf = zeros(nB, ny)
vd = zeros(1, ny)
vc = zeros(nB, ny)
policy = zeros(nB, ny)
q = ones(nB, ny) .* (1 / (1 + r))
defprob = zeros(nB, ny)
return (β = β, γ = γ, r = r, ρ = ρ, η = η, θ = θ, ny = ny,
nB = nB, ygrid = ygrid, ydefgrid = ydefgrid,
Bgrid = Bgrid, Π = Π, vf = vf, vd = vd, vc = vc,
policy = policy, q = q, defprob = defprob)
return (;beta, gamma, r, rho, eta, theta, ny,
nB, ygrid, ydefgrid, Bgrid, Pi, vf, vd, vc,
policy, q, defprob)
end
u(ae, c) = c^(1 - ae.γ) / (1 - ae.γ)
u(ae, c) = c^(1 - ae.gamma) / (1 - ae.gamma)
function one_step_update!(ae,
EV,
EVd,
EVc)
# unpack stuff
@unpack β, γ, r, ρ, η, θ, ny, nB = ae
@unpack ygrid, ydefgrid, Bgrid, Π, vf, vd, vc, policy, q, defprob = ae
(;beta, gamma, r, rho, eta, theta, ny, nB) = ae
(;ygrid, ydefgrid, Bgrid, Pi, vf, vd, vc, policy, q, defprob) = ae
zero_ind = searchsortedfirst(Bgrid, 0.)
for iy in 1:ny
y = ae.ygrid[iy]
ydef = ae.ydefgrid[iy]
# value of being in default with income y
defval = u(ae, ydef) + β * (θ * EVc[zero_ind, iy] + (1-θ) * EVd[1, iy])
defval = u(ae, ydef) + beta * (theta * EVc[zero_ind, iy] + (1-theta) * EVd[1, iy])
ae.vd[1, iy] = defval
for ib in 1:nB
Expand All @@ -371,7 +370,7 @@ function one_step_update!(ae,
pol_ind = 0
for ib_next=1:nB
c = max(y - ae.q[ib_next, iy]*Bgrid[ib_next] + B, 1e-14)
m = u(ae, c) + β * EV[ib_next, iy]
m = u(ae, c) + beta * EV[ib_next, iy]
if m > current_max
current_max = m
Expand All @@ -390,24 +389,24 @@ end
function compute_prices!(ae)
# unpack parameters
@unpack β, γ, r, ρ, η, θ, ny, nB = ae
(;beta, gamma, r, rho, eta, theta, ny, nB) = ae
# create default values with a matching size
vd_compat = repeat(ae.vd, nB)
default_states = vd_compat .> ae.vc
# update default probabilities and prices
copyto!(ae.defprob, default_states * ae.Π')
copyto!(ae.defprob, default_states * ae.Pi')
copyto!(ae.q, (1 .- ae.defprob) / (1 + r))
return
end
function vfi!(ae; tol = 1e-8, maxit = 10000)
# unpack stuff
@unpack β, γ, r, ρ, η, θ, ny, nB = ae
@unpack ygrid, ydefgrid, Bgrid, Π, vf, vd, vc, policy, q, defprob = ae
Πt = Π'
(;beta, gamma, r, rho, eta, theta, ny, nB) = ae
(;ygrid, ydefgrid, Bgrid, Pi, vf, vd, vc, policy, q, defprob) = ae
Pit = Pi'
# Iteration stuff
it = 0
Expand All @@ -420,11 +419,11 @@ function vfi!(ae; tol = 1e-8, maxit = 10000)
it += 1
# compute expectations for this iterations
# (we need Π' because of order value function dimensions)
# (we need Pi' because of order value function dimensions)
copyto!(V_upd, ae.vf)
EV = ae.vf * Πt
EVd = ae.vd * Πt
EVc = ae.vc * Πt
EV = ae.vf * Pit
EVd = ae.vd * Pit
EVc = ae.vc * Pit
# update value function
one_step_update!(ae, EV, EVd, EVc)
Expand Down Expand Up @@ -452,7 +451,7 @@ function QuantEcon.simulate(ae,
B_init_ind = searchsortedfirst(ae.Bgrid, B_init)
# create a QE MarkovChain
mc = MarkovChain(ae.Π)
mc = MarkovChain(ae.Pi)
y_sim_indices = simulate(mc, capT + 1; init = y_init_ind)
# allocate and fill output
Expand Down Expand Up @@ -495,8 +494,8 @@ function QuantEcon.simulate(ae,
y_sim_val[t] = ae.ydefgrid[y_sim_indices[t]]
q_sim_val[t] = ae.q[zero_index, y_sim_indices[t]]
# with probability θ exit default status
default_status[t + 1] = rand() ≥ ae.θ
# with probability theta exit default status
default_status[t + 1] = rand() >= ae.theta
end
end
Expand Down Expand Up @@ -591,12 +590,12 @@ using DataFrames, Plots
Compute the value function, policy and equilibrium prices

```{code-cell} julia
ae = ArellanoEconomy(β = .953, # time discount rate
γ = 2., # risk aversion
ae = ArellanoEconomy(beta = .953, # time discount rate
gamma = 2., # risk aversion
r = 0.017, # international interest rate
ρ = .945, # persistence in output
η = 0.025, # st dev of output shock
θ = 0.282, # prob of regaining access
rho = .945, # persistence in output
eta = 0.025, # st dev of output shock
theta = 0.282, # prob of regaining access
ny = 21, # number of points in y grid
nB = 251) # number of points in B grid
Expand Down Expand Up @@ -627,7 +626,7 @@ q_low = zeros(0)
q_high = zeros(0)
for i in 1:ae.nB
b = ae.Bgrid[i]
if -0.35 ≤ b ≤ 0 # to match fig 3 of Arellano
if -0.35 <= b <= 0 # to match fig 3 of Arellano
push!(x, b)
push!(q_low, ae.q[i, iy_low])
push!(q_high, ae.q[i, iy_high])
Expand Down
22 changes: 11 additions & 11 deletions lectures/multi_agent_models/harrison_kreps.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -290,12 +290,12 @@ Here's a function that can be used to compute these values
using LinearAlgebra
function price_single_beliefs(transition, dividend_payoff;
β=.75)
beta=.75)
# First compute inverse piece
imbq_inv = inv(I - β * transition)
imbq_inv = inv(I - beta * transition)
# Next compute prices
prices = β * ((imbq_inv * transition) * dividend_payoff)
prices = beta * ((imbq_inv * transition) * dividend_payoff)
return prices
end
Expand Down Expand Up @@ -416,19 +416,19 @@ Here's code to solve for $\bar p$, $\hat p_a$ and $\hat p_b$ using the iterative
```{code-cell} julia
function price_optimistic_beliefs(transitions,
dividend_payoff;
β=.75, max_iter=50000,
beta=.75, max_iter=50000,
tol=1e-16)
# We will guess an initial price vector of [0, 0]
p_new = [0,0]
p_old = [10.0,10.0]
# We know this is a contraction mapping, so we can iterate to conv
for i ∈ 1:max_iter
for i in 1:max_iter
p_old = p_new
temp = [maximum((q * p_old) + (q * dividend_payoff))
for q in transitions]
p_new = β * temp
p_new = beta * temp
# If we succed in converging, break out of for loop
if maximum(sqrt, ((p_new - p_old).^2)) < 1e-12
Expand All @@ -437,7 +437,7 @@ function price_optimistic_beliefs(transitions,
end
temp=[minimum((q * p_old) + (q * dividend_payoff)) for q in transitions]
ptwiddle = β * temp
ptwiddle = beta * temp
phat_a = [p_new[1], ptwiddle[2]]
phat_b = [ptwiddle[1], p_new[2]]
Expand Down Expand Up @@ -486,17 +486,17 @@ Here's code to solve for $\check p$ using iteration
```{code-cell} julia
function price_pessimistic_beliefs(transitions,
dividend_payoff;
β=.75, max_iter=50000,
beta=.75, max_iter=50000,
tol=1e-16)
# We will guess an initial price vector of [0, 0]
p_new = [0, 0]
p_old = [10.0, 10.0]
# We know this is a contraction mapping, so we can iterate to conv
for i ∈ 1:max_iter
for i in 1:max_iter
p_old = p_new
temp=[minimum((q * p_old) + (q* dividend_payoff)) for q in transitions]
p_new = β * temp
p_new = beta * temp
# If we succed in converging, break out of for loop
if maximum(sqrt, ((p_new - p_old).^2)) < 1e-12
Expand Down Expand Up @@ -594,7 +594,7 @@ heterogeneous beliefs.
opt_beliefs = price_optimistic_beliefs([qa, qb], dividendreturn)
labels = ["p_optimistic", "p_hat_a", "p_hat_b"]
for (p, label) ∈ zip(opt_beliefs, labels)
for (p, label) in zip(opt_beliefs, labels)
println(label)
println(repeat("=", 20))
s0, s1 = round.(p, digits = 2)
Expand Down
Loading

0 comments on commit befd388

Please sign in to comment.