/
index.qmd
236 lines (185 loc) · 5.81 KB
/
index.qmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
---
title: "Scientific Communication with Quarto: Economic Models"
date: 2023-06-10
last-updated: last-modified
author: Kazuharu Yanagimoto
categories: [Julia, Observable JS, Quarto]
image: "img/anim.gif"
twitter-card:
image: "img/anim.gif"
execute:
warning: false
message: false
format:
html:
code-fold: true
code-tools: true
---
## Quarto Is not Just for Static Reports
Notebook coding (Jupyter, Rmarkdown, and Quarto) has gotten popular in the last decade.
It allows us to tell a story with code and data, and it makes scientific communication more transparent and reproducible.
I use Quarto to share my current results with my supervisor and coauthors,
especially when I am in the playground stage of the project.
While it is naturally useful in that it is an easy tool for creating reports,
exploiting its HTML feature, we can enrich communication even further.
Here, I introduce my use case of GIF and Observable JS for economic models with Julia.
Imagine that you are working on a very simple life-cycle model:
For $t = 1, \dots, T$, households solve
$$
V(t, e, x) = \max_{c, x'} \frac{c^{1 - \sigma}}{1 - \sigma} + \beta \mathbb{E}V(t + 1, e', x')
$$
$$
\begin{aligned}
c + x' &\le (1 + r)x + ew \\
\text{Pr}(e' | e) &= \Gamma(e) \\
c, x' &\ge 0
\end{aligned}
$$
See the following code for the parameters and the solution algorithm.
```{julia}
using QuantEcon
using Plots
using LaTeXStrings
using CSV
using DataFrames
@kwdef struct Model
# Grid for x
nₓ::Int = 30
x̲::Float64 = 0.1
x̄::Float64 = 4.0
x_grid::Vector{Float64} = range(start = x̲, stop = x̄, length = nₓ)
# Grid for e: parameters for Tauchen
nₑ::Int = 10
σ_ε::Float64 = 0.02058
λ_ε::Float64 = 0.99
m::Float64 = 1.5
mc::MarkovChain = tauchen(nₑ, λ_ε, σ_ε, 0.0, m)
e_grid::Vector{Float64} = exp.(mc.state_values)
P::Matrix{Float64} = mc.p
# Utility function
σ::Float64 = 2.0
β::Float64 = 0.97
T::Int = 10
# Prices
r::Float64 = 0.07
w::Float64 = 5.0
end
u(c, m::Model) = isone(m.σ) ? log(c) : c^(1 - m.σ) / (1 - m.σ);
function solve(m::Model)
(; T, nₓ, nₑ, r, w, β, P, x_grid, e_grid) = m
V = zeros(nₓ, nₑ, T)
for t = T:-1:1, iₓ= 1:nₓ, iₑ = 1:nₑ
utility = -Inf
for iₓ′ = 1:nₓ
expected = (t == T) ? 0.0 :
sum(P[iₑ, iₑ′] * V[iₓ′, iₑ′, t+1] for iₑ′ = 1:nₑ)
c = (1 + r) * x_grid[iₓ] + e_grid[iₑ] * w - x_grid[iₓ′]
if c > 0
utility = max(u(c, m) + β * expected, utility)
end
end
V[iₓ, iₑ, t] = utility
end
return V
end
m = Model()
V = solve(m);
```
A traditional way to show the value function is to plot for selected $e$ and $t$.
```{julia}
ps = []
y̲, ȳ = minimum(V) * 1.1, maximum(V) + 0.1
for (i, t) ∈ enumerate([1, 4, 7, 10])
p = plot(m.x_grid, V[:, 1, t],
xlabel = "x",
ylims = (y̲, ȳ),
label = L"e_1",
legend = :bottomright,
title = "t = $t")
plot!(m.x_grid, V[:, 5, t], label = L"e_5")
plot!(m.x_grid, V[:, 10, t], label = L"e_{10}")
push!(ps, p)
end
plot(ps...)
```
This is not bad, however, we can add more information by using GIF.
## GIF
While this is not so famous, we can create a GIF animation just by adding `@animate` macro in Julia.
```{julia}
#| label: thumbnail
#| code-fold: show
anim = @animate for t = 1:10
plot(m.x_grid, V[:, 1, t],
xlabel = "x",
ylims = (y̲, ȳ),
label = L"e_1",
legend = :bottomright,
title = "t = $t")
plot!(m.x_grid, V[:, 5, t], label = L"e_5")
plot!(m.x_grid, V[:, 10, t], label = L"e_{10}")
end
gif(anim, "img/anim.gif", fps = 1)
```
See? This is very intuitive and informative!
But if you are interested in the changes with multiple parameters, you may want to try Observable JS.
## Observable JS
[Observable](https://observablehq.com/) is a JavaScript notebook that allows us to create interactive visualizations.
It is free and open-source, and we can embed it in our Quarto document.
Since the Observable is JavaScript-based, we can put an interactive (dynamic) visualization in **static** HTML documents!
For example, you are now interested in the changes in the value function with $r$ and $w$.
To implement it in Observable, we first need to create a CSV file^[You can use CSV, TSV, JSON, Arrow (uncompressed), and SQLite as data input. See [Data Sources](https://quarto.org/docs/computations/ojs.html#data-sources)] that contains the results of the simulation.
```{julia}
#| code-fold: show
function solve_partial(r, w)
m = Model(r = r, w = w)
V = solve(m)
return [
(x = x, ie = iₑ, t = t, V = V[iₓ, iₑ, t], r = m.r, w = m.w)
for (iₓ, x) ∈ enumerate(m.x_grid)
for (iₑ, e) ∈ enumerate(m.e_grid) if iₑ ∈ [1, 5, 10]
for t ∈ 1:m.T
]
end
res = Iterators.flatten([solve_partial(r, w) for r ∈ 0.01:0.02:0.15 for w ∈ 2.0:0.5:10.0])
df = DataFrame(res)
d = Dict(1 => "e₁", 5 => "e₅", 10 => "e₁₀")
df.lbl = [d[i] for i in df.ie]
CSV.write("data.csv", df);
```
This CSV file is then loaded in Observable and filtered by the parameters.^[While the parameters `r`, `w`, and `t` are defined in the next cell, the order of the Observable JS cells do not change the results.]
```{ojs}
//| code-fold: show
data = FileAttachment("data.csv").csv()
filtered = data.filter(function(sim) {
return sim.r == r && sim.w == w && sim.t == t
})
```
Adding sliders and a plot, we have a cool interactive visualization!
```{ojs}
//| panel: input
viewof t = Inputs.range(
[1, 10],
{value: 1, step: 1, label: "t"}
)
viewof r = Inputs.range(
[0.01, 0.15],
{value: 0.07, step: 0.02, label: "r"}
)
viewof w = Inputs.range(
[2.0, 10.0],
{value: 5.0, step: 0.5, label: "w"}
)
```
```{ojs}
Plot.plot({
marginLeft: 50,
height: 400,
color: {domain: ["e₁", "e₅", "e₁₀"], legend: true},
x: {domain: [0.0, 4.0]},
y: {domain: [-5.5, 0.0]},
marks: [
Plot.lineY(filtered, {x: "x", y: "V", stroke: "lbl"}),
]
})
```
Happy Quarto life 🥂!