The source files for all examples can be found in /examples.
Example 4: Pareto surface
This example kicks up the complexity a couple of notches. We will introduce a new optimisation estimator, NearOptimalCentering optimiser.
using PortfolioOptimisers, PrettyTables
# Format for pretty tables.
tsfmt = (v, i, j) -> begin
if j == 1
return Date(v)
else
return v
end
end;
resfmt = (v, i, j) -> begin
if j == 1
return v
else
return isa(v, Number) ? "$(round(v*100, digits=3)) %" : v
end
end;1. ReturnsResult data
We will use the same data as the previous example.
using CSV, TimeSeries, DataFrames
X = TimeArray(CSV.File(joinpath(@__DIR__, "SP500.csv.gz")); timestamp = :Date)[(end - 252):end]
pretty_table(X[(end - 5):end]; formatters = tsfmt)
# Compute the returns
rd = prices_to_returns(X)ReturnsResult
nx | 20-element Vector{String}
X | 252×20 Matrix{Float64}
nf | nothing
F | nothing
ts | 252-element Vector{Dates.Date}
iv | nothing
ivpa | nothing
2. Preparing solvers for pareto surface
The pareto surface is a generalisation of the efficient frontier, in fact, we can even think of hypersurfaces if we provide more parameters, but that would be difficult to visualise, so we will stick to a 2D surface in 3D space.
We'll provide a vector of solvers beacause the optimisation type we'll be using is more complex, and will contain various constraints.
using Clarabel
slv = [Solver(; name = :clarabel1, solver = Clarabel.Optimizer,
settings = Dict("verbose" => false),
check_sol = (; allow_local = true, allow_almost = true)),
Solver(; name = :clarabel2, solver = Clarabel.Optimizer,
settings = Dict("verbose" => false, "max_step_fraction" => 0.95),
check_sol = (; allow_local = true, allow_almost = true)),
Solver(; name = :clarabel3, solver = Clarabel.Optimizer,
settings = Dict("verbose" => false, "max_step_fraction" => 0.9),
check_sol = (; allow_local = true, allow_almost = true)),
Solver(; name = :clarabel4, solver = Clarabel.Optimizer,
settings = Dict("verbose" => false, "max_step_fraction" => 0.85),
check_sol = (; allow_local = true, allow_almost = true)),
Solver(; name = :clarabel5, solver = Clarabel.Optimizer,
settings = Dict("verbose" => false, "max_step_fraction" => 0.8),
check_sol = (; allow_local = true, allow_almost = true)),
Solver(; name = :clarabel6, solver = Clarabel.Optimizer,
settings = Dict("verbose" => false, "max_step_fraction" => 0.75),
check_sol = (; allow_local = true, allow_almost = true)),
Solver(; name = :clarabel7, solver = Clarabel.Optimizer,
settings = Dict("verbose" => false, "max_step_fraction" => 0.70),
check_sol = (; allow_local = true, allow_almost = true))];3. High order prior statistics
We will once again precompute the prior statistics because otherwise they'd have to be recomputed a few times.
We will be using high order risk measures, so we need to compute high order moments, we can do this with a HighOrderPriorEstimator, which needs a prior estimator that computes low order moments. Since we are only using a year of data, we will denoise our positive definite matrices by eliminating the eigenvalues corresponding to random noise. Denoising the non-positive definite matrix for the data we're using creates a negative square root, so we will not denoise it.
Note how many options this estimator contains.
de = Denoise(; alg = SpectralDenoise(;))
mp = DefaultMatrixProcessing(; denoise = de)
pe = HighOrderPriorEstimator(;
# Prior estimator for low order moments
pe = EmpiricalPrior(;
ce = PortfolioOptimisersCovariance(;
mp = mp)),
# Estimator for cokurtosis
kte = Cokurtosis(; mp = mp),
# Estimator for coskewness
ske = Coskewness())HighOrderPriorEstimator
pe | EmpiricalPrior
| ce | PortfolioOptimisersCovariance
| | ce | Covariance
| | | me | SimpleExpectedReturns
| | | | w | nothing
| | | ce | GeneralWeightedCovariance
| | | | ce | StatsBase.SimpleCovariance: StatsBase.SimpleCovariance(true)
| | | | w | nothing
| | | alg | Full()
| | mp | DefaultMatrixProcessing
| | | pdm | Posdef
| | | | alg | UnionAll: NearestCorrelationMatrix.Newton
| | | denoise | Denoise
| | | | alg | SpectralDenoise()
| | | | args | Tuple{}: ()
| | | | kwargs | @NamedTuple{}: NamedTuple()
| | | | kernel | typeof(AverageShiftedHistograms.Kernels.gaussian): AverageShiftedHistograms.Kernels.gaussian
| | | | m | Int64: 10
| | | | n | Int64: 1000
| | | detone | nothing
| | | alg | nothing
| me | SimpleExpectedReturns
| | w | nothing
| horizon | nothing
kte | Cokurtosis
| me | SimpleExpectedReturns
| | w | nothing
| mp | DefaultMatrixProcessing
| | pdm | Posdef
| | | alg | UnionAll: NearestCorrelationMatrix.Newton
| | denoise | Denoise
| | | alg | SpectralDenoise()
| | | args | Tuple{}: ()
| | | kwargs | @NamedTuple{}: NamedTuple()
| | | kernel | typeof(AverageShiftedHistograms.Kernels.gaussian): AverageShiftedHistograms.Kernels.gaussian
| | | m | Int64: 10
| | | n | Int64: 1000
| | detone | nothing
| | alg | nothing
| alg | Full()
ske | Coskewness
| me | SimpleExpectedReturns
| | w | nothing
| mp | NonPositiveDefiniteMatrixProcessing
| | denoise | nothing
| | detone | nothing
| | alg | nothing
| alg | Full()
Lets compute the prior statistics.
pr = prior(pe, rd)HighOrderPrior
pr | LowOrderPrior
| X | 252×20 Matrix{Float64}
| mu | 20-element Vector{Float64}
| sigma | 20×20 Matrix{Float64}
| chol | nothing
| w | nothing
| ens | nothing
| kld | nothing
| ow | nothing
| rr | nothing
| f_mu | nothing
| f_sigma | nothing
| f_w | nothing
kt | 400×400 Matrix{Float64}
L2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
S2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
sk | 20×400 Matrix{Float64}
V | 20×20 Matrix{Float64}
skmp | NonPositiveDefiniteMatrixProcessing
| denoise | nothing
| detone | nothing
| alg | nothing
In order to generate a pareto surface/hyper-surface, we need more dimensions than we've previously explored. We can do this by adding more risk measure sweeps (and taking their product) to generate a mesh. PortfolioOptimisers does this internally and generally, but we will limit ourselves to two risk measures. This will generate a 2D surface which we can visualise in 3D.
We will use the square root NegativeSkewness and SquareRootKurtosis.
r1 = NegativeSkewness()
r2 = SquareRootKurtosis()SquareRootKurtosis
settings | RiskMeasureSettings
| scale | Float64: 1.0
| ub | nothing
| rke | Bool: true
w | nothing
mu | nothing
kt | nothing
N | nothing
alg | Full()
4. Near optimal centering pareto surface
First we need to get the bounds of our pareto surface. We can do this in many different ways, the simplest are:
- Minimise the risk using both risk measures simultaneously subject to optional constraints.
- Maximise the return, utility or ratio subject to optional constraints.
We will simply maximise the risk-return ratio for both risk measures on their own with no added constraints. This will not give a complete surface, but it will give us a reasonable range of values.
The NearOptimalCentering estimator will not return the portfolio which satisfies the traditional MeanRisk constraints, but rather a portfolio which is at the centre of an analytical region (neighbourhood) around the optimal solution. The region is parametrised by binning the efficient frontier, we will use the automatic bins here, but it is possible to define them manually.
# Risk-free rate of 4.2/100/252
rf = 4.2 / 100 / 252
opt = JuMPOptimiser(; pe = pr, slv = slv)
obj = MaximumRatio(; rf = rf)
opt1 = NearOptimalCentering(; r = r1, obj = obj, opt = opt)
opt2 = NearOptimalCentering(; r = r2, obj = obj, opt = opt)NearOptimalCentering
opt | JuMPOptimiser
| pe | HighOrderPrior
| | pr | LowOrderPrior
| | | X | 252×20 Matrix{Float64}
| | | mu | 20-element Vector{Float64}
| | | sigma | 20×20 Matrix{Float64}
| | | chol | nothing
| | | w | nothing
| | | ens | nothing
| | | kld | nothing
| | | ow | nothing
| | | rr | nothing
| | | f_mu | nothing
| | | f_sigma | nothing
| | | f_w | nothing
| | kt | 400×400 Matrix{Float64}
| | L2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
| | S2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
| | sk | 20×400 Matrix{Float64}
| | V | 20×20 Matrix{Float64}
| | skmp | NonPositiveDefiniteMatrixProcessing
| | | denoise | nothing
| | | detone | nothing
| | | alg | nothing
| slv | 7-element Vector{Solver{Symbol, UnionAll, T3, @NamedTuple{allow_local::Bool, allow_almost::Bool}, Bool} where T3}
| wb | WeightBounds
| | lb | Float64: 0.0
| | ub | Float64: 1.0
| bgt | Float64: 1.0
| sbgt | nothing
| lt | nothing
| st | nothing
| lcs | nothing
| lcm | nothing
| cent | nothing
| gcard | nothing
| sgcard | nothing
| smtx | nothing
| sgmtx | nothing
| slt | nothing
| sst | nothing
| sglt | nothing
| sgst | nothing
| sets | nothing
| nplg | nothing
| cplg | nothing
| tn | nothing
| te | nothing
| fees | nothing
| ret | ArithmeticReturn
| | ucs | nothing
| | lb | nothing
| sce | SumScalariser: SumScalariser()
| ccnt | nothing
| cobj | nothing
| sc | Int64: 1
| so | Int64: 1
| card | nothing
| scard | nothing
| nea | nothing
| l1 | nothing
| l2 | nothing
| ss | nothing
| strict | Bool: false
r | SquareRootKurtosis
| settings | RiskMeasureSettings
| | scale | Float64: 1.0
| | ub | nothing
| | rke | Bool: true
| w | nothing
| mu | nothing
| kt | nothing
| N | nothing
| alg | Full()
obj | MaximumRatio
| rf | Float64: 0.0001666666666666667
| ohf | nothing
bins | nothing
w_min | nothing
w_min_ini | nothing
w_opt | nothing
w_opt_ini | nothing
w_max | nothing
w_max_ini | nothing
ucs_flag | Bool: true
alg | UnconstrainedNearOptimalCentering()
Note the number of options in the estimator. In particular the alg property. Which in this case means the NearOptimalCentering alg will not have any external constraints applied to it.
Lets optimise the portfolios.
res1 = optimise!(opt1)
res2 = optimise!(opt2)NearOptimalCenteringOptimisation
oe | DataType: NearOptimalCentering
pa | ProcessedJuMPOptimiserAttributes
| pr | HighOrderPrior
| | pr | LowOrderPrior
| | | X | 252×20 Matrix{Float64}
| | | mu | 20-element Vector{Float64}
| | | sigma | 20×20 Matrix{Float64}
| | | chol | nothing
| | | w | nothing
| | | ens | nothing
| | | kld | nothing
| | | ow | nothing
| | | rr | nothing
| | | f_mu | nothing
| | | f_sigma | nothing
| | | f_w | nothing
| | kt | 400×400 Matrix{Float64}
| | L2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
| | S2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
| | sk | 20×400 Matrix{Float64}
| | V | 20×20 Matrix{Float64}
| | skmp | NonPositiveDefiniteMatrixProcessing
| | | denoise | nothing
| | | detone | nothing
| | | alg | nothing
| wb | WeightBounds
| | lb | 20-element StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}
| | ub | 20-element StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}
| lt | nothing
| st | nothing
| lcs | nothing
| cent | nothing
| gcard | nothing
| sgcard | nothing
| smtx | nothing
| sgmtx | nothing
| slt | nothing
| sst | nothing
| sglt | nothing
| sgst | nothing
| nplg | nothing
| cplg | nothing
| tn | nothing
| fees | nothing
| ret | ArithmeticReturn
| | ucs | nothing
| | lb | nothing
w_min_retcode | OptimisationSuccess
| res | Dict{Any, Any}: Dict{Any, Any}()
w_opt_retcode | OptimisationSuccess
| res | Dict{Any, Any}: Dict{Any, Any}()
w_max_retcode | OptimisationSuccess
| res | Dict{Any, Any}: Dict{Any, Any}()
noc_retcode | OptimisationSuccess
| res | Dict{Any, Any}: Dict{Any, Any}(:clarabel1 => Dict{Symbol, Any}(:settings => Dict{String, Bool}("verbose" => 0), :err => solution_summary(; result = 1, verbose = false)
| ├ solver_name : Clarabel
| ├ Termination
| │ ├ termination_status : SLOW_PROGRESS
| │ ├ result_count : 1
| │ └ raw_status : INSUFFICIENT_PROGRESS
| ├ Solution (result = 1)
| │ ├ primal_status : OTHER_RESULT_STATUS
| │ ├ dual_status : OTHER_RESULT_STATUS
| │ ├ objective_value : 9.29400e+01
| │ └ dual_objective_value : 9.29050e+01
| └ Work counters
| ├ solve_time (sec) : 2.61538e-01
| └ barrier_iterations : 35))
retcode | OptimisationSuccess
| res | nothing
sol | PortfolioOptimisers.JuMPOptimisationSolution
| w | 20-element Vector{Float64}
model | A JuMP Model
| ├ solver: Clarabel
| ├ objective_sense: MIN_SENSE
| │ └ objective_function_type: JuMP.AffExpr
| ├ num_variables: 273
| ├ num_constraints: 47
| │ ├ JuMP.AffExpr in MOI.EqualTo{Float64}: 1
| │ ├ Vector{JuMP.AffExpr} in MOI.Nonnegatives: 1
| │ ├ Vector{JuMP.AffExpr} in MOI.Nonpositives: 1
| │ ├ Vector{JuMP.AffExpr} in MOI.SecondOrderCone: 1
| │ ├ Vector{JuMP.AffExpr} in MOI.ExponentialCone: 42
| │ └ Vector{JuMP.AffExpr} in MOI.PositiveSemidefiniteConeSquare: 1
| └ Names registered in the model
| └ :M_PSD, :W, :W_M, :bgt, :ckurt_soc_1, :clog_delta_w, :clog_ret, :clog_risk, :clog_w, :k, :log_delta_w, :log_ret, :log_risk, :log_w, :lw, :obj_expr, :ret, :risk, :risk_vec, :sc, :so, :sqrt_kurtosis_risk_1, :w, :w_lb, :w_ub, :zkurt_1
In order to allow for multiple risk measures in optimisations, certain measures can take different parameters. In this case, NegativeSkewness and SquareRootKurtosis take the moment matrices, which are used to compute the risk measures. We can use the factory function to create a new risk measure with the same parameters as the original, but with the moment matrices from the prior. Other risk measures require a solver, and this function is also used in those cases.
r1 = factory(r1, pr)
r2 = factory(r2, pr)SquareRootKurtosis
settings | RiskMeasureSettings
| scale | Float64: 1.0
| ub | nothing
| rke | Bool: true
w | nothing
mu | 20-element Vector{Float64}
kt | 400×400 Matrix{Float64}
N | nothing
alg | Full()
Lets compute the risk bounds for the pareto surface. We need to compute four risks because we have two risk measures and two optimisations. This will let us pick the lower and upper bounds for each risk measure, as we explore the pareto surface from one optimisation to the other.
sk_rk1 = expected_risk(r1, res1.w, pr.X);
kt_rk1 = expected_risk(r2, res1.w, pr.X);
sk_rk2 = expected_risk(r1, res2.w, pr.X);
kt_rk2 = expected_risk(r2, res2.w, pr.X);We will now create new risk measures bounded by these values. We will also use factories from the get-go. The optimisation procedure prioritises the parameters in the risk measures over the ones in the prior. This lets users provide the same risk measure with different parameters in the same optimisation. We will use two ranges of 5. The total number of points in the pareto surface will be the product of the points of each range.
Since we don't know which sk_rk1 or sk_r2, kt_rk1 or kt_rk2 is bigger or smaller, we need to use min, max.
r1 = factory(NegativeSkewness(;
settings = RiskMeasureSettings(;
# Risk upper bounds go from the minimum to maximum risk given the optimisations.
ub = range(;
start = min(sk_rk1,
sk_rk2),
stop = max(sk_rk1,
sk_rk2),
length = 5))), pr);
r2 = factory(SquareRootKurtosis(;
settings = RiskMeasureSettings(;
ub = range(;
start = min(kt_rk1,
kt_rk2),
stop = max(kt_rk1,
kt_rk2),
length = 5))), pr);Now we only need to maximise the return given both risk measures. Internally, the optimisation will generate the mesh as a product of the ranges in the order in which the risk measures were provided. This also works with the MeanRisk estimatro, in fact, NearOptimalCentering uses it internally.
Since we are using an unconstrained NearOptimalCentering, the risk bound constraints will not be satisfied by the solution. If we wish to satisfy them, we can provide alg = ConstrainedNearOptimalCentering(), but would also make the optimisations harder, which may cause them to fail.
opt3 = NearOptimalCentering(; r = [r1, r2], obj = MaximumReturn(), opt = opt)NearOptimalCentering
opt | JuMPOptimiser
| pe | HighOrderPrior
| | pr | LowOrderPrior
| | | X | 252×20 Matrix{Float64}
| | | mu | 20-element Vector{Float64}
| | | sigma | 20×20 Matrix{Float64}
| | | chol | nothing
| | | w | nothing
| | | ens | nothing
| | | kld | nothing
| | | ow | nothing
| | | rr | nothing
| | | f_mu | nothing
| | | f_sigma | nothing
| | | f_w | nothing
| | kt | 400×400 Matrix{Float64}
| | L2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
| | S2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
| | sk | 20×400 Matrix{Float64}
| | V | 20×20 Matrix{Float64}
| | skmp | NonPositiveDefiniteMatrixProcessing
| | | denoise | nothing
| | | detone | nothing
| | | alg | nothing
| slv | 7-element Vector{Solver{Symbol, UnionAll, T3, @NamedTuple{allow_local::Bool, allow_almost::Bool}, Bool} where T3}
| wb | WeightBounds
| | lb | Float64: 0.0
| | ub | Float64: 1.0
| bgt | Float64: 1.0
| sbgt | nothing
| lt | nothing
| st | nothing
| lcs | nothing
| lcm | nothing
| cent | nothing
| gcard | nothing
| sgcard | nothing
| smtx | nothing
| sgmtx | nothing
| slt | nothing
| sst | nothing
| sglt | nothing
| sgst | nothing
| sets | nothing
| nplg | nothing
| cplg | nothing
| tn | nothing
| te | nothing
| fees | nothing
| ret | ArithmeticReturn
| | ucs | nothing
| | lb | nothing
| sce | SumScalariser: SumScalariser()
| ccnt | nothing
| cobj | nothing
| sc | Int64: 1
| so | Int64: 1
| card | nothing
| scard | nothing
| nea | nothing
| l1 | nothing
| l2 | nothing
| ss | nothing
| strict | Bool: false
r | Vector{PortfolioOptimisers.RiskMeasure}: PortfolioOptimisers.RiskMeasure[NegativeSkewness
settings | RiskMeasureSettings
| scale | Float64: 1.0
| ub | StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}: 0.0014407874342906597:0.00010182732908744173:0.0018480967506404266
| rke | Bool: true
mp | NonPositiveDefiniteMatrixProcessing
| denoise | nothing
| detone | nothing
| alg | nothing
sk | 20×400 Matrix{Float64}
V | 20×20 Matrix{Float64}
alg | SqrtRiskExpr()
, SquareRootKurtosis
settings | RiskMeasureSettings
| scale | Float64: 1.0
| ub | StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}: 0.00023763521983569873:7.758896959421372e-6:0.0002686708076733842
| rke | Bool: true
w | nothing
mu | 20-element Vector{Float64}
kt | 400×400 Matrix{Float64}
N | nothing
alg | Full()
]
obj | MaximumReturn()
bins | nothing
w_min | nothing
w_min_ini | nothing
w_opt | nothing
w_opt_ini | nothing
w_max | nothing
w_max_ini | nothing
ucs_flag | Bool: true
alg | UnconstrainedNearOptimalCentering()
See how r is a vector of risk measures with populated properties. We can now optimise the porftolios.
res3 = optimise!(opt3)NearOptimalCenteringOptimisation
oe | DataType: NearOptimalCentering
pa | ProcessedJuMPOptimiserAttributes
| pr | HighOrderPrior
| | pr | LowOrderPrior
| | | X | 252×20 Matrix{Float64}
| | | mu | 20-element Vector{Float64}
| | | sigma | 20×20 Matrix{Float64}
| | | chol | nothing
| | | w | nothing
| | | ens | nothing
| | | kld | nothing
| | | ow | nothing
| | | rr | nothing
| | | f_mu | nothing
| | | f_sigma | nothing
| | | f_w | nothing
| | kt | 400×400 Matrix{Float64}
| | L2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
| | S2 | 210×400 SparseArrays.SparseMatrixCSC{Int64, Int64}
| | sk | 20×400 Matrix{Float64}
| | V | 20×20 Matrix{Float64}
| | skmp | NonPositiveDefiniteMatrixProcessing
| | | denoise | nothing
| | | detone | nothing
| | | alg | nothing
| wb | WeightBounds
| | lb | 20-element StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}
| | ub | 20-element StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}
| lt | nothing
| st | nothing
| lcs | nothing
| cent | nothing
| gcard | nothing
| sgcard | nothing
| smtx | nothing
| sgmtx | nothing
| slt | nothing
| sst | nothing
| sglt | nothing
| sgst | nothing
| nplg | nothing
| cplg | nothing
| tn | nothing
| fees | nothing
| ret | ArithmeticReturn
| | ucs | nothing
| | lb | nothing
w_min_retcode | OptimisationSuccess
| res | Dict{Any, Any}: Dict{Any, Any}()
w_opt_retcode | 25-element Vector{PortfolioOptimisers.OptimisationReturnCode}
w_max_retcode | OptimisationSuccess
| res | Dict{Any, Any}: Dict{Any, Any}()
noc_retcode | 25-element Vector{PortfolioOptimisers.OptimisationReturnCode}
retcode | OptimisationSuccess
| res | nothing
sol | 25-element Vector{PortfolioOptimisers.JuMPOptimisationSolution}
model | A JuMP Model
| ├ solver: Clarabel
| ├ objective_sense: MIN_SENSE
| │ └ objective_function_type: JuMP.AffExpr
| ├ num_variables: 274
| ├ num_constraints: 48
| │ ├ JuMP.AffExpr in MOI.EqualTo{Float64}: 1
| │ ├ Vector{JuMP.AffExpr} in MOI.Nonnegatives: 1
| │ ├ Vector{JuMP.AffExpr} in MOI.Nonpositives: 1
| │ ├ Vector{JuMP.AffExpr} in MOI.SecondOrderCone: 2
| │ ├ Vector{JuMP.AffExpr} in MOI.ExponentialCone: 42
| │ └ Vector{JuMP.AffExpr} in MOI.PositiveSemidefiniteConeSquare: 1
| └ Names registered in the model
| └ :M_PSD, :W, :W_M, :bgt, :ckurt_soc_2, :clog_delta_w, :clog_ret, :clog_risk, :clog_w, :cnskew_soc_1, :k, :log_delta_w, :log_ret, :log_risk, :log_w, :lw, :nskew_risk_1, :obj_expr, :ret, :risk, :risk_vec, :sc, :so, :sqrt_kurtosis_risk_2, :w, :w_lb, :w_ub, :zkurt_2
As expected, there are 5 × 5 = 25 solutions. Thankfully there are no warnings about failed optimisations, so there is no need to check the solutions.
The NearOptimalCentering estimator contains various return codes because it may need to compute some MeanRisk optimisations, it has a retcode which summarises whether all other optimisations succeeded. We can check this to make sure it was a success.
isa(res3.retcode, OptimisationSuccess)true5. Visualising the pareto surface
Lets view how the weights evolve along the pareto surface.
using StatsPlots, GraphRecipes
plot_stacked_area_composition(res3.w, rd.nx)
Now we can view the parteo surface. For the z-axis and colourbar, we will use the conditional drawdown at risk to return ratio.
plot_measures(res3.w, pr; x = r1, y = r2,
z = RatioRiskMeasure(; rk = ConditionalDrawdownatRisk(),
rt = ArithmeticReturn(), rf = rf),
c = RatioRiskMeasure(; rk = ConditionalDrawdownatRisk(),
rt = ArithmeticReturn(), rf = rf),
title = "Pareto Surface", xlabel = "Sqrt NSke", ylabel = "Sqrt Kt",
zlabel = "CDaR/Return")
We can view it in 2D as well.
gr()
plot_measures(res3.w, pr; x = r1, y = r2,
c = RatioRiskMeasure(; rk = ConditionalDrawdownatRisk(),
rt = ArithmeticReturn(), rf = rf),
title = "Pareto Front", xlabel = "Sqrt NSke", ylabel = "Sqrt Kt",
colorbar_title = "\n\nCDaR/Return", right_margin = 8Plots.mm)
This page was generated using Literate.jl.