Nearest correlation

This example illustrates the sensitivity analysis of the nearest correlation problem studied in [H02].

Higham, Nicholas J. Computing the nearest correlation matrix—a problem from finance. IMA journal of Numerical Analysis 22.3 (2002): 329-343.

using DiffOpt, JuMP, SCS, LinearAlgebra
solver = SCS.Optimizer

function proj(A, dH = Diagonal(ones(size(A, 1))), H = ones(size(A)))
    n = LinearAlgebra.checksquare(A)
    model = Model(() -> DiffOpt.diff_optimizer(solver))
    @variable(model, X[1:n, 1:n] in PSDCone())
    @constraint(model, [i in 1:n], X[i, i] == 1)
    @objective(model, Min, sum((H .* (X - A)) .^ 2))
    MOI.set(
        model,
        DiffOpt.ForwardObjectiveFunction(),
        sum((dH .* (X - A)) .^ 2),
    )
    optimize!(model)
    DiffOpt.forward_differentiate!(model)
    dX = MOI.get.(model, DiffOpt.ForwardVariablePrimal(), X)
    return value.(X), dX
end

proj (generic function with 3 methods)

Example from [H02, p. 334-335]:

A = LinearAlgebra.Tridiagonal(ones(2), ones(3), ones(2))

3×3 LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}:
 1.0  1.0   ⋅ 
 1.0  1.0  1.0
  ⋅   1.0  1.0

The projection is computed as follows:

X, dX = proj(A)

------------------------------------------------------------------
	       SCS v3.3.1 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 6, constraints m: 9
cones: 	  z: primal zero / dual free vars: 3
	  s: psd vars: 6, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 5
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 9, nnz(P): 6
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 9.98e-01  3.10e+00  2.08e+01 -1.22e+01  1.00e-01  1.98e-04
    50| 1.65e-04  8.63e-06  3.40e-04 -6.72e+00  1.00e-01  5.40e-04
------------------------------------------------------------------
status:  solved
timings: total: 5.41e-04s = setup: 9.56e-05s + solve: 4.45e-04s
	 lin-sys: 1.81e-05s, cones: 3.04e-04s, accel: 4.09e-06s
------------------------------------------------------------------
objective = -6.721638
------------------------------------------------------------------

The projection of A is:

X

3×3 Matrix{Float64}:
 1.0       0.760853  0.157204
 0.760853  1.0       0.760853
 0.157204  0.760853  1.0

The derivative of the projection with respect to a uniform increase of the weights of the diagonal entries is:

dX

3×3 Matrix{Float64}:
 -1.03268e-7  -0.0648562   0.0240023
 -0.0648562   -3.0611e-7  -0.0648563
  0.0240023   -0.0648563   2.49923e-7

Example from [H02, Section 4, p. 340]:

A = LinearAlgebra.Tridiagonal(-ones(3), 2ones(4), -ones(3))

4×4 LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}:
  2.0  -1.0    ⋅     ⋅ 
 -1.0   2.0  -1.0    ⋅ 
   ⋅   -1.0   2.0  -1.0
   ⋅     ⋅   -1.0   2.0

The projection is computed as follows:

X, dX = proj(A)

------------------------------------------------------------------
	       SCS v3.3.1 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 10, constraints m: 14
cones: 	  z: primal zero / dual free vars: 4
	  s: psd vars: 10, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 5
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 14, nnz(P): 10
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 1.01e+00  3.32e+00  2.85e+01 -3.31e+01  1.00e-01  1.51e-04
    50| 1.87e-05  1.19e-06  5.04e-05 -1.74e+01  1.00e-01  5.63e-04
------------------------------------------------------------------
status:  solved
timings: total: 5.64e-04s = setup: 7.86e-05s + solve: 4.86e-04s
	 lin-sys: 2.39e-05s, cones: 3.82e-04s, accel: 4.68e-06s
------------------------------------------------------------------
objective = -17.447231
------------------------------------------------------------------

The projection of A is:

X

4×4 Matrix{Float64}:
  1.0       -0.808422   0.191579   0.106771
 -0.808422   1.0       -0.656252   0.191579
  0.191579  -0.656252   1.0       -0.808422
  0.106771   0.191579  -0.808422   1.0

The derivative of the projection with respect to a uniform increase of the weights of the diagonal entries is:

dX

4×4 Matrix{Float64}:
 -9.36251e-8  -0.174219   -0.0684906   -0.031452
 -0.174219    -3.3989e-7  -0.250254    -0.0684904
 -0.0684906   -0.250254    3.06644e-7  -0.174219
 -0.031452    -0.0684904  -0.174219     2.89576e-7

This page was generated using Literate.jl.