# Modeling dynamic economic systems using ModelingToolkit (MTK)

I am seeking help in specifying a dynamic economic model in MTK. Although I have followed the evolution of MTK since its creation, I have difficulties with “function” as a function can be the whole model, a function of parameters or a function of variables. My aim is to simulate the dynamic model written in MTK to generate synthetic timeseries and then take it across to DataDrivenDiffEq for model discovery.

The nonlinear model consists of a system of three first order ODEs:

\begin{aligned} \frac{\dot{ω}}{ω} & = \Big(\frac{η}{1+η}\Big) \Big[Φ(λ) −a_l\Big] \ \ (1) \\ \frac{\dot{λ}}{λ} & = \Big[κ(π) C \Big(\frac{1 - ω}{b}\Big)^\frac{1}{η} −δ−a_l−β - \frac{1}{η(1-ω)} \frac{\dot{ω}}{ω}\Big] \ \ (2) \\ \frac{\dot{d_f}}{d_f} & = \Big[r − g (κ(π), ω)\Big] + \frac{κ(π) −(1 − ω)}{d_f} \ \ (3) \\ \end{aligned}

with four behavioural/variable functions π, Φ(λ) , κ(π) and g(κ(π), ω) of the form:

\begin{aligned} π & : = 1 - ω - r d_f \\ Φ(λ) & = Φ₁ + Φ₂e^{(Φ₃λ)} \\ κ(π) & = κ₁ + κ₂e^{(κ₃π)} \\ g(κ(π), ω) & = κ(π) C \Big(\frac{1 - ω}{b}\Big)^\frac{1}{η} -δ - \frac{\dot{ω}}{(1 - ω)η} \\ \end{aligned}

the parameters are: t η κ C b δ β a Φ₁ Φ₂ Φ₃ κ₁ κ₂ κ₃ r

and the variables are: ω(t), λ(t), d_f (t)

This nonlinear formulation is already relatively simple as it is specified in logarithmic derivative terms. Otherwise it is even more complicated and nonlinear. Help in defining and registering the various types of functions would be much appreciated.

1 Like

I don’t understand the question. What’s your current code?

using LinearAlgebra
using ModelingToolkit
using OrdinaryDiffEq
using Plots
using Random
using Symbolics: scalarize
using Latexify

Random.seed!(1111) # For the noise

@parameters t η κ C b δ β a π Φ₁ Φ₂ Φ₃ κ₁ κ₂ κ₃ r
@variables ω(t) λ(t) d_f(t)
D = Differential(t)

function π(ω, r, d_f) # Accounting identity
1 - ω - r*d_f
end
@register π

function Φ(λ, Φ₁, Φ₂, Φ₃) # Phillip’s curve
Φ₁ + Φ₂exp(Φ₃λ)
end
@register Φ(λ)

function κ(π, κ₁, κ₂, κ₃) # Investment function
κ₁ + κ₂exp(κ₃π)
end
@register κ(π)

function g(κ, ω) # Growth rate of the economy
κC((1-ω)/b)^(1/η) - δ - D(ω)/((1 - ω)*η)
end
@ register g(κ, ω)

# Define a nonlinear system

sfc = [D(ω) ~ ω*(η/(1+η))(Φ(λ) - a),
D(λ) ~ λ
(κ(π)C((1-ω)/b)^(1/η) - δ -β -a -(1/(η*(1-ω)))D(ω)/ω),
D(d_f) ~ d_f
(r - g(κ, ω)) + κ(π) - (1 - ω)]

latexify(sfc)

That seems like it would just work if you didn’t register the functions? You want to register as little as possible so it can define symbolic derivatives by tracing.

For now, I get the right math for the model with this:

@variables t ω(t) λ(t) Φ(t) κ(t) π(t) d_f(t) g(t) # independent and dependent variables
@parameters Φ₁ Φ₂ Φ₃ η a_l κ₁ κ₂ κ₃ δ β r C b # parameters
D = Differential(t) # define an operator for the differentiation w.r.t. time

sfc = [ Φ ~ Φ₁ + Φ₂exp(Φ₃λ),
κ ~ κ₁ + κ₂exp(κ₃π),
π ~ (1 - ω - rd_f),
g ~ κ
C*((1-ω)/b)^(1/η) - δ - D(ω)/((1 - ω)η),
D(ω) ~ ω
((η/(1+η))(Φ - a_l)),
D(λ) ~ λ
((κC((1-ω)/b)^(1/η) - δ -β -a_l -(1/(η*(1-ω)))D(ω)/ω)),
D(d_f) ~ d_f
(r - g) + κ - (1 - ω)]

latexify(sfc)

The next step is to see how “NonlinearSystem” and “NonlinearProblem” will accept this sfc model. If it solves as it is I do not need functions at all and this reads as textbook presentation.

That’s not a nonlinear system. It’s an ODESystem representing a DAE, which you can solve with DAEProblem's constructor.

As this is economics, isn’t this a system derived from an optimization problem with state, control and co-state variable? And the derivatives are with respect to time? And you have initial values and a terminal/limit condition? And do you know the nature of this system, e.g. saddle-point stable?

I am just lost in the process of setting up the DAEProblem. The docs and test of the Robertson DAE problem is too cryptic for me. Any other DAEProblem examples using ModelingToolkit anywhere?

This is where I stand:
@variables t ω(t) λ(t) Φ(t) κ(t) π(t) d_f(t) g(t) # independent and dependent variables
@parameters Φ₁ Φ₂ Φ₃ η a_l κ₁ κ₂ κ₃ δ β r C b # parameters
D = Differential(t) # define an operator for the differentiation w.r.t. time

sfc = [ res1 = Φ₁ + Φ₂exp(Φ₃λ) - Φ,
res2 = κ₁ + κ₂exp(κ₃π) - κ,
res3 = (1 - ω - rd_f) - π,
res4 = κ
C*((1-ω)/b)^(1/η) - δ - D(ω)/((1 - ω)η) - g,
res5 = ω
((η/(1+η))(Φ - a_l)) - D(ω),
res6 = λ
((κC((1-ω)/b)^(1/η) - δ -β -a_l -(1/(η*(1-ω)))D(ω)/ω)) - D(λ),
res7 = d_f
(r - g) + κ - (1 - ω) - D(d_f)]

latexify(sfc)

and then
@named keen = NonlinearSystem(sfc, t, [ω, λ, d_f, Φ, κ, π, g], [Φ₁, Φ₂, Φ₃, η, a_l, κ₁, κ₂, κ₃, δ, β, r, C, b])

u0 = [ω => 0.75,
λ => 0.75,
d_f => 0.80]
du0 = [ ]

param = [Φ₁ => -0.01, Φ₂ => , Φ₃ => ,
η => , a_l => 2.0 ,
κ₁ => , κ₂ => , κ₃ => ,
δ => 1, β => 1, r => 4.0, C => 0.333, b => 0.135]

prob = DAEProblem(keen,du0,u0,(0.0,50.0), param)
sol = solve(prob,IDA())

Not yet. For now it is just a stock-flow consistent model of three nonlinear ODEs that can exhibit chaotic behaviour for some values of parameters and initial values. Will get fancier when this basic setup works.

Nice to see a Steve Keen fan here. Are you considering to develop an SFC package in Julia using MLT or just exploring?

See the DAE benchmarks for examples, like Chemical Akzo Nobel Differential-Algebraic Equation (DAE) Work-Precision Diagrams

Great!

Let me give the solution to the MTK-DAE problem for the interested readers. This could be simplified further with a D ln x (t) operator as mentioned in Logarithmic differentiation in MTK · Issue #1344 · SciML/ModelingToolkit.jl · GitHub

using LinearAlgebra
using ModelingToolkit
using DifferentialEquations: solve
using NonlinearSolve
using OrdinaryDiffEq
using Plots: plot
using Random
using Symbolics: scalarize
using Latexify

ModelingToolkit.@parameters begin
Φ₁=0.04/(1. - 0.04^2)
Φ₂=0.04^3/(1. - 0.04^2)
η=500
a₁=0.02
κ₁=0.
κ₂=0.05
κ₃=1.75
δ=0.05
β=0.01
r=0.03
C=0.33
b=0.135
end

@variables begin
t
ω(t) = 0.75
λ(t) = 0.90
d₁(t) = 0.5
π(t) = 0.14
Y(t) = 100.0
Φ(t) # = 0.95
κ(t) # = 0.05
g(t) # = 0.05
end

D = Differential(t)

Φ = -Φ₁ + Φ₂/((1. - λ)^2)
κ = κ₁ + κ₂exp(κ₃π)
g = κC((1 - ω)/b)^(1/η) - δ - ω/((1 - ω)(1 + η))(Φ - a₁)

sfc = [
D(ω) ~ ω*((η/(1 + η))(Φ - a₁))
D(λ) ~ λ
((κC((1 - ω)/b)^(1/η) - δ -β -a₁ -(1/((1 + η)(1 - ω)))(Φ - a₁)))
D(d₁) ~ d₁*(r - g) + κ - (1 - ω)
D(Y) ~ Yg
0. ~ π - 1 + ω + r
d₁
]

latexify(sfc)

ModelingToolkit.@named keen = ModelingToolkit.ODESystem(sfc)

tspan = (0.0, 200.0)
mmprob = ODEProblem(keen, , tspan)
sol = solve(mmprob, Rodas4(),abstol=1/10^14,reltol=1/10^14);

du = mmprob.f(mmprob.u0,mmprob.p,0.0)
du0 = D.(states(keen)) .=> du
daeprob = DAEProblem(keen,du0,,tspan)
ref_sol = solve(daeprob,IDA(),abstol=1/10^14,reltol=1/10^14);

#probs = [mmprob,daeprob]
#refs = [ref_sol,ref_sol];

The idea is to replicate in Julia using MTK the SFC models in Giraud and Grasselli models (https://ms.mcmaster.ca/~grasselli/GiraudGrasselli2021.pdf) and in Bastidas et al (https://www.parisschoolofeconomics.eu/docs/fabre-adrien/minskyan-classical-growth-cycles--bastidas,-fabre,-mcisaac-2018.pdf).

My roadmap is to show that continuous time dynamical systems in economics that have been around since the 1970s can be discovered and identified from the dataset using DataDrivenDiffEq. This is a two step process: 1) SFC model producing synthetic data => data discovery of SFC model and 2) use real data for the same variables to discover which model is compatible with it.

1 Like

Sounds nice.
Gaël Giraud is a friend’s thesis supervisor, let me know if you want to be in touch.