Skip to content

Commit

Permalink
Merge pull request #6 from thorek1/models_wo_past_or_shocks
Browse files Browse the repository at this point in the history 
Models wo past or shocks
  • Loading branch information
@thorek1
thorek1 authored Dec 2, 2022
2 parents 9db33cc + 71cde2d commit 8018762
Show file tree
Hide file tree
Showing 4 changed files with 69 additions and 15 deletions.
44 changes: 30 additions & 14 deletions src/MacroModelling.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1365,9 +1365,14 @@ function calculate_jacobian(parameters::Vector{<: Number}, 𝓂::ℳ)
par = ComponentVector( vcat(parameters,calibrated_parameters),Axis(vcat(𝓂.parameters,𝓂.calibration_equations_parameters)))
SS = ComponentVector(non_stochastic_steady_state, Axis(sort(union(𝓂.exo_present,𝓂.var))))

SS_past = SS[[indexin(sort([𝓂.var_past; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_past,𝓂.exo_past))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_past))...]
SS_present = SS[[indexin(sort([𝓂.var_present; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_present,𝓂.exo_present))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_present))...]
SS_future = SS[[indexin(sort([𝓂.var_future; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_future,𝓂.exo_future))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_future))...]
past_idx = [indexin(sort([𝓂.var_past; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_past,𝓂.exo_past))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_past = length(past_idx) > 0 ? SS[past_idx] : zeros(0) #; zeros(length(𝓂.exo_past))...]

present_idx = [indexin(sort([𝓂.var_present; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_present,𝓂.exo_present))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_present = length(present_idx) > 0 ? SS[present_idx] : zeros(0)#; zeros(length(𝓂.exo_present))...]

future_idx = [indexin(sort([𝓂.var_future; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_future,𝓂.exo_future))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_future = length(future_idx) > 0 ? SS[future_idx] : zeros(0)#; zeros(length(𝓂.exo_future))...]

shocks_ss = zeros(length(𝓂.exo))

Expand All @@ -1383,10 +1388,15 @@ function calculate_hessian(parameters::Vector{<: Number}, 𝓂::ℳ)

par = ComponentVector( vcat(parameters,calibrated_parameters),Axis(vcat(𝓂.parameters,𝓂.calibration_equations_parameters)))
SS = ComponentVector(non_stochastic_steady_state, Axis(𝓂.var))

SS_past = SS[[indexin(sort([𝓂.var_past; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_past,𝓂.exo_past))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_past))...]
SS_present = SS[[indexin(sort([𝓂.var_present; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_present,𝓂.exo_present))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_present))...]
SS_future = SS[[indexin(sort([𝓂.var_future; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_future,𝓂.exo_future))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_future))...]

past_idx = [indexin(sort([𝓂.var_past; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_past,𝓂.exo_past))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_past = length(past_idx) > 0 ? SS[past_idx] : zeros(0) #; zeros(length(𝓂.exo_past))...]

present_idx = [indexin(sort([𝓂.var_present; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_present,𝓂.exo_present))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_present = length(present_idx) > 0 ? SS[present_idx] : zeros(0)#; zeros(length(𝓂.exo_present))...]

future_idx = [indexin(sort([𝓂.var_future; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_future,𝓂.exo_future))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_future = length(future_idx) > 0 ? SS[future_idx] : zeros(0)#; zeros(length(𝓂.exo_future))...]

shocks_ss = zeros(length(𝓂.exo))

Expand All @@ -1404,10 +1414,15 @@ function calculate_third_order_derivatives(parameters::Vector{<: Number}, 𝓂::

par = ComponentVector( vcat(parameters,calibrated_parameters),Axis(vcat(𝓂.parameters,𝓂.calibration_equations_parameters)))
SS = ComponentVector(non_stochastic_steady_state, Axis(𝓂.var))

SS_past = SS[[indexin(sort([𝓂.var_past; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_past,𝓂.exo_past))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_past))...]
SS_present = SS[[indexin(sort([𝓂.var_present; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_present,𝓂.exo_present))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_present))...]
SS_future = SS[[indexin(sort([𝓂.var_future; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_future,𝓂.exo_future))]), sort(union(𝓂.var,𝓂.exo_present)))...]]#; zeros(length(𝓂.exo_future))...]

past_idx = [indexin(sort([𝓂.var_past; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_past,𝓂.exo_past))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_past = length(past_idx) > 0 ? SS[past_idx] : zeros(0) #; zeros(length(𝓂.exo_past))...]

present_idx = [indexin(sort([𝓂.var_present; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_present,𝓂.exo_present))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_present = length(present_idx) > 0 ? SS[present_idx] : zeros(0)#; zeros(length(𝓂.exo_present))...]

future_idx = [indexin(sort([𝓂.var_future; map(x -> Symbol(replace(string(x), r"ᴸ⁽⁻[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾|ᴸ⁽[⁰¹²³⁴⁵⁶⁷⁸⁹]+⁾" => "")), union(𝓂.aux_future,𝓂.exo_future))]), sort(union(𝓂.var,𝓂.exo_present)))...]
SS_future = length(future_idx) > 0 ? SS[future_idx] : zeros(0)#; zeros(length(𝓂.exo_future))...]

shocks_ss = zeros(length(𝓂.exo))

Expand Down Expand Up @@ -1868,10 +1883,11 @@ function irf(state_update::Function, initial_state::Vector{Float64}, T::timings;
return KeyedArray(Y[var_idx,:,:]; Variables = T.var[var_idx], Period = 1:periods, Shock = [:simulate])
elseif shocks == :none
Y = zeros(T.nVars,periods,1)
Y[:,1,1] = state_update(initial_state,[0.0])
shck = T.nExo == 0 ? Vector{Float64}(undef, 0) : [0.0]
Y[:,1,1] = state_update(initial_state,shck)

for t in 1:periods-1
Y[:,t+1,1] = state_update(Y[:,t,1],[0.0])
Y[:,t+1,1] = state_update(Y[:,t,1],shck)
end

return KeyedArray(Y[var_idx,:,:]; Variables = T.var[var_idx], Period = 1:periods, Shock = [:none])
Expand Down Expand Up @@ -1993,7 +2009,7 @@ function parse_shocks_input_to_index(shocks::Symbol_input, T::timings)
if shocks == :all
shock_idx = 1:T.nExo
elseif shocks == :none
shock_idx = 1:T.nExo
shock_idx = 1
elseif shocks == :simulate
shock_idx = 1
elseif shocks isa Matrix{Symbol}
Expand Down
8 changes: 8 additions & 0 deletions src/get_functions.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -62,6 +62,8 @@ function get_irf(𝓂::ℳ,
sol_mat = calculate_first_order_solution(jacc; T = 𝓂.timings)

state_update = function(state::Vector, shock::Vector) sol_mat * [state[𝓂.timings.past_not_future_and_mixed_idx]; shock] end

shocks = 𝓂.timings.nExo == 0 ? :none : shocks

shock_idx = parse_shocks_input_to_index(shocks,𝓂.timings)

Expand Down Expand Up @@ -181,6 +183,12 @@ function get_irf(𝓂::ℳ;
var_idx = parse_variables_input_to_index(variables, 𝓂.timings)
end

shocks = 𝓂.timings.nExo == 0 ? :none : shocks

if shocks == :none && generalised_irf
@error "Cannot compute generalised IRFs for model without shocks."
end

if generalised_irf
girfs = girf(state_update,
𝓂.timings;
Expand Down
4 changes: 3 additions & 1 deletion src/plotting.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -76,6 +76,8 @@ function plot(𝓂::ℳ;
end

init_state = initial_state == [0.0] ? zeros(𝓂.timings.nVars) : initial_state - collect(get_non_stochastic_steady_state_internal(𝓂))

shocks = 𝓂.timings.nExo == 0 ? :none : shocks

shock_idx = parse_shocks_input_to_index(shocks,𝓂.timings)

Expand Down Expand Up @@ -206,7 +208,6 @@ function plot(𝓂::ℳ;
# elseif length(pp) > 1
if length(pp) > 0

shock_string = ": " * string(𝓂.timings.exo[shock_idx[shock]])

if shocks == :simulate
shock_string = ": simulate all"
Expand All @@ -215,6 +216,7 @@ function plot(𝓂::ℳ;
shock_string = ""
shock_name = "no_shock"
else
shock_string = ": " * string(𝓂.timings.exo[shock_idx[shock]])
shock_name = string(𝓂.timings.exo[shock_idx[shock]])
end

Expand Down
28 changes: 28 additions & 0 deletions test/runtests.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,34 @@ using Random
include("test_standalone_function.jl")
end

@testset "Model without shocks" begin
@model m begin
K[0] = (1 - δ) * K[-1] + I[0]
Z[0] = (1 - ρ) * μ + ρ * Z[-1]
I[1] = ((ρ + δ - Z[0])/(1 - δ)) + ((1 + ρ)/(1 - δ)) * I[0]
end

@parameters m begin
ρ = 0.05
δ = 0.10
μ = .17
σ = .2
end

m_ss = get_steady_state(m)
@test isapprox(m_ss(:,:Steady_state),[1/7.5,1/.75,.17],rtol = eps(Float32))

m_sol = get_solution(m)
@test isapprox(m_sol(:,:K),[1/.75,.9,.04975124378109454],rtol = eps(Float32))

init = m_ss(:,:Steady_state) |> collect
init[2] *= 1.5
get_irf(m, initial_state = init, shocks = :none)

plot_irf(m, initial_state = init, shocks = :none)
@test true
end

@testset "Distribution functions, general and SS" begin

@model RBC_CME begin
Expand Down

0 comments on commit 8018762

Please sign in to comment.