using CairoMakie
using ColorSchemes
using Distributions
using LabelledArrays
using Random
Replicating the AAA equtity generator
modeling
scenario-generator
This notebook replicates the model and parameters for the real world equity generator described in this AAA 2005 reference paper.
Stochastic Log Volatility Model
Note that the @.
and other broadcasting (.
symbol) allows us to operate on multiple funds at once.
function v(v_prior,params,Zₜ)
(;σ_v, σ_m,σ_p,σ⃰,ϕ,τ) = params
v_m = log.(σ_m)
v_p = log.(σ_p)
v⃰ = log.(σ⃰)
# vol are the odd values in the random array
ṽ = @. min(v_p, (1 - ϕ) * v_prior + ϕ * log(τ) ) + σ_v * Zₜ[[1,3,5,7]]
v = @. max(v_m, min(v⃰,ṽ))
return v
end
function scenario(params,Z;months=1200)
(;σ_v,σ_0, ρ,A,B,C) = params
n_funds = size(params,2)
#initilize/pre-allocate
Zₜ = rand(Z)
v_t = log.(σ_0)
σ_t = zeros(n_funds)
μ_t = zeros(n_funds)
log_returns = map(1:months) do t
Zₜ = rand!(Z,Zₜ)
v_t .= v(v_t,params,Zₜ)
σ_t .= exp.(v_t)
@. μ_t = A + B * σ_t + C * (σ_t)^2
# return are the even values in the random array
log_return = @. μ_t / 12 + σ_t / sqrt(12) * Zₜ[[2,4,6,8]]
end
# convert vector of vector to matrix
reduce(hcat,log_returns)
end
scenario (generic function with 1 method)
Model Parameters
# use a labelled array for easy reference of the parameters
params = @LArray [
0.12515 0.14506 0.16341 0.20201 # τ
0.35229 0.41676 0.3632 0.35277 # ϕ
0.32645 0.32634 0.35789 0.34302 # σ_v
-0.2488 -0.1572 -0.2756 -0.2843 # ρ
0.055 0.055 0.055 0.055 # A
0.56 0.466 0.67 0.715 # B
-0.9 -0.9 -0.95 -1.0 # C
0.1476 0.1688 0.2049 0.2496 # σ_0
0.0305 0.0354 0.0403 0.0492 # σ_m
0.3 0.3 0.4 0.55 # σ_p
0.7988 0.4519 0.9463 1.1387 # σ⃰
] (
# define the regions each label refers to
τ = (1,:),
ϕ = (2,:),
σ_v = (3,:),
ρ = (4,:),
A = (5,:),
B = (6,:),
C = (7,:),
σ_0 = (8,:),
σ_m = (9,:),
σ_p = (10,:),
σ⃰ = (11,:)
)
11×4 LArray{Float64, 2, Matrix{Float64}, (τ = (1, Colon()), ϕ = (2, Colon()), σ_v = (3, Colon()), ρ = (4, Colon()), A = (5, Colon()), B = (6, Colon()), C = (7, Colon()), σ_0 = (8, Colon()), σ_m = (9, Colon()), σ_p = (10, Colon()), σ⃰ = (11, Colon()))}:
:τ => 0.12515 … :τ => 0.20201
:ϕ => 0.35229 :ϕ => 0.35277
:σ_v => 0.32645 :σ_v => 0.34302
:ρ => -0.2488 :ρ => -0.2843
:A => 0.055 :A => 0.055
:B => 0.56 … :B => 0.715
:C => -0.9 :C => -1.0
:σ_0 => 0.1476 :σ_0 => 0.2496
:σ_m => 0.0305 :σ_m => 0.0492
:σ_p => 0.3 :σ_p => 0.55
:σ⃰ => 0.7988 … :σ⃰ => 1.1387
The Multivariate normal and covariance matrix
# 11 columns because it's got the bond returns in it
cov_matrix = [
1.000 -0.249 0.318 -0.082 0.625 -0.169 0.309 -0.183 0.023 0.075 0.080;
-0.249 1.000 -0.046 0.630 -0.123 0.829 -0.136 0.665 -0.120 0.192 0.393;
0.318 -0.046 1.000 -0.157 0.259 -0.050 0.236 -0.074 -0.066 0.034 0.044;
-0.082 0.630 -0.157 1.000 -0.063 0.515 -0.098 0.558 -0.105 0.130 0.234;
0.625 -0.123 0.259 -0.063 1.000 -0.276 0.377 -0.180 0.034 0.028 0.054;
-0.169 0.829 -0.050 0.515 -0.276 1.000 -0.142 0.649 -0.106 0.067 0.267;
0.309 -0.136 0.236 -0.098 0.377 -0.142 1.000 -0.284 0.026 0.006 0.045;
-0.183 0.665 -0.074 0.558 -0.180 0.649 -0.284 1.000 0.034 -0.091 -0.002;
0.023 -0.120 -0.066 -0.105 0.034 -0.106 0.026 0.034 1.000 0.047 -0.028;
0.075 0.192 0.034 0.130 0.028 0.067 0.006 -0.091 0.047 1.000 0.697;
0.080 0.393 0.044 0.234 0.054 0.267 0.045 -0.002 -0.028 0.697 1.000;
]
Z = MvNormal(
zeros(11), #means for return and volatility
cov_matrix # covariance matrix
# full covariance matrix in AAA Excel workook on Parameters tab
)
FullNormal(
dim: 11
μ: [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
Σ: [1.0 -0.249 … 0.075 0.08; -0.249 1.0 … 0.192 0.393; … ; 0.075 0.192 … 1.0 0.697; 0.08 0.393 … 0.697 1.0]
)
Scenarios and validation
A single scenario
x = scenario(params,Z;months=1200)
4×1200 Matrix{Float64}:
-0.0124354 0.00912116 0.0214326 … 0.0641472 -0.0224623 -0.000456449
0.0225229 0.0423898 0.0161604 0.0540821 -0.0368879 0.0440957
-0.0655741 0.022554 0.0250733 0.0955047 -0.0303337 0.00775239
-0.144574 0.0262148 -0.00133166 0.0932583 -0.144408 0.0400758
Validation of summary statistics
The summary statistics expected (per paper Table 8):
μ ≈ [0.0060, 0.0062, 0.0063, 0.0065]
σ ≈ [0.0436, 0.0492, 0.0590, 0.0724]
These computed values match very closely:
# generate 1000 scenarios
scens = [scenario(params,Z) for _ in 1:1000];
let
# compute summary statistics
μ = vec(mean(mean(x,dims=2) for x in scens))
σ = vec(mean(std(x,dims=2) for x in scens))
(;μ,σ)
end
(μ = [0.006045006563819091, 0.0062199437086345784, 0.006278168390979234, 0.006531685069702403], σ = [0.043517357354043225, 0.04913425078838221, 0.0590681180411829, 0.07212209999154577])
Plotting some scenarios
let
f = Figure()
n = 25
colors = ColorSchemes.Johnson
ax = Axis(f[1,1],yscale=log10,ylabel="index value",
title="$n realizations of 4 correlated equity funds per AAA ESG")
for s in scens[1:n]
for i in 1:4
lines!(ax,cumprod(exp.(s[i,:])), color=(colors[i],0.3),label="fund $i")
end
end
axislegend(ax,unique=true,position=:lt)
f
end
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227