Replicating the AAA equtity generator

using CairoMakie
using ColorSchemes
using Distributions
using LabelledArrays
using Random

include(joinpath(@__DIR__, "..", "..", "assets", "themes", "cassette_futurism.jl"))
set_theme!(cassette_futurism_theme())

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
    ṽ =  @. min(v_p, (1 - ϕ) * v_prior + ϕ * log(τ) ) + σ_v * Zₜ[[1,3,5,7]]
    
    v = @. max(v_m, min(v⃰,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 LabelledArrays.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.0226621  0.0943664  0.0234509   …  0.050131   0.0492229  -0.0219075
 -0.0555582  0.0829307  0.00521493     0.0339099  0.0828788   0.0431517
  0.0316437  0.132638   0.0219535      0.0543943  0.0816699  -0.0176026
 -0.0140118  0.0759929  0.0408952      0.0308992  0.0359973   0.00172504

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.006026485565774503, 0.006136467570645955, 0.006263016508235973, 0.0064386917541197995], σ = [0.04366248954214508, 0.049159561611857854, 0.05911064823506265, 0.07223410658240541])

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