# 1 - Setting up
using Plots
α = 0.3
δ = 0.9
φ = α/(1-α*δ)
ψ = 1/(1-δ)*(log(1-α*δ)+(α*δ/(1-α*δ))*log(α*δ))
# 1.a
N = 101
k_min = 0
k_max = 10
step = (k_max-k_min)/(N-1)
K = [i for i in k_min:step:k_max]
# 1.b
v(k) = ψ + φ*log(k)
V = [v(k) for k in K]
# You can get the same result by doing this
V = v.(K) # It works better when you only have one argument.
# 1.c
plot(K,V) # Change the value of N above to create different plots.
# The lower the value of N, the worse this approximation is.
# Part 2
# 2.a
# If you coded the Tauchen method in math camp you can use it now. The
# Tauchen method will, in effect, set Py = 0.1 for all values except the
# extremes, which will have probability 0.05 each.
Y = [y for y in 0:0.1:1]
Py = ones(11,1).*0.1
Py[1] = 0.05
Py[11] = 0.05
# 2.b.
K = [k for k in 0:0.1:9]
k_vals = repeat(K,1,11)
y_vals = repeat(Y,1,91)'
new_vals = k_vals + y_vals
# 2.c
V_vals = V[Int.(round.((new_vals .+0.1)*10))] # Use round if Int alone doesn't work.
# 2.d
new_V = V_vals * Py
# 2.e
utility_increase = abs.((new_V - V[1:91])./V[1:91]) # The percentage increase isn't quite defined given that we have a negative utility, but this works
plot(K,utility_increase)
# 2.f
# We used k less than 9 because otherwise we'd fall out of the grid. One
# potential solution is to interpolate values outside of the grid (see
# problem 6). Another one is to assume that V(k>10) = V(10)
function interpolate(x,X,F)
N = length(F)
if x < minimum(X) # Not relevant in this case since it'll be -Inf
return F[1] + (x-X[1])/(X[2]-X[1])*(F[2]-F[1]);
end
if x >= maximum(X)
return F[N] + (x-X[N])/(X[N]-X[N-1]) * (F[N]-F[N-1]);
end
for i in 1:N-1
if x == X[i]
return F[i];
elseif x >= X[i] && x < X[i+1]
return F[i] + (x-X[i])/(X[i+1]-X[i]) * (F[i+1]-F[i]);
end
end
end
K = [i for i in k_min:step:k_max]
V = [v(k) for k in K]
Y = [y for y in 0:0.01:1]
Py = ones(101,1).*0.1
Py[1] = 0.05
Py[101] = 0.05
# 2.b.
k_vals = repeat(K,1,101)
y_vals = repeat(Y,1,101)'
new_vals = k_vals + y_vals
# 2.c
V_vals = zeros(101,101)
for i in 1:101
for j in 1:101
V_vals[i,j] = interpolate(new_vals[i,j],K,V)
end
end
# 2.d
new_V = V_vals * Py