Next Article in Journal
Fixed Circle and Fixed Disc Results for New Types of Θc-Contractive Mappings in Metric Spaces
Next Article in Special Issue
Utilizing Macro Fiber Composite to Control Rotating Blade Vibrations
Previous Article in Journal
Multi-Scale Study on Mechanical Property and Strength of New Green Sand (Poly Lactic Acid) as Replacement of Fine Aggregate in Concrete Mix
Previous Article in Special Issue
Some New Oscillation Results for Fourth-Order Neutral Differential Equations with Delay Argument
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Computing Nearest Correlation Matrix via Low-Rank ODE’s Based Technique

by
Mutti-Ur Rehman
1,2,
Jehad Alzabut
3 and
Kamaleldin Abodayeh
3,*
1
Department of Mathematics, Sukkur IBA University, Sukkur 65200, Pakistan
2
Department of Chemical Engineering, Massachusetts Institute of Technology, Cambridge, MA 02139, USA
3
Department of Mathematics and General Sciences, Prince Sultan University, Riyadh 11586, Saudi Arabia
*
Author to whom correspondence should be addressed.
Symmetry 2020, 12(11), 1824; https://doi.org/10.3390/sym12111824
Submission received: 16 October 2020 / Revised: 31 October 2020 / Accepted: 2 November 2020 / Published: 4 November 2020

Abstract

:
For n-dimensional real-valued matrix A, the computation of nearest correlation matrix; that is, a symmetric, positive semi-definite, unit diagonal and off-diagonal entries between 1 and 1 is a problem that arises in the finance industry where the correlations exist between the stocks. The proposed methodology presented in this article computes the admissible perturbation matrix and a perturbation level to shift the negative spectrum of perturbed matrix to become non-negative or strictly positive. The solution to optimization problems constructs a gradient system of ordinary differential equations that turn over the desired perturbation matrix. Numerical testing provides enough evidence for the shifting of the negative spectrum and the computation of nearest correlation matrix.

1. Introduction

The correlation matrices A i R n , n are real, symmetric A i t = A i , positive semi-definite λ j ( A i ) 0 , j , with λ j being the spectrum of A i and having a unit diagonal, d i a g ( A i ) = 1 . These correlation matrices appear when one is interested in computing correlations between pairs of random variables. For example, in the finance industry, the correlation between stocks measured over a fixed amount of time, and the missing data can compromise the correlations and hence will direct to a non-positive semi-definite matrix. Similarly, a practitioner can explore the effect on a portfolio that assigns the correlations between assets differently from a certain amount of historical values, which can easily destroy the semi-definiteness of the matrix (or matrices); for more details, refer to [1,2,3].
The approximation of correlation matrices has led much intention and interest in finding a nearest correlation matrix X that minimizes the matrix Frobenius norm of A X in matrix nearness problem
m i n { A X F : X = X t , X 0 , d i a g ( X ) = 1 } ,
where for X , Y being as the symmetric matrices and the quantity X Y is positive semi-definite matrix and X F = t r a c e ( X t X ) as R n , n represents Hilbert space with X , Y = t r a c e ( X t Y ) . The constraints in above minimization problem are closed and convex sets so that the matrix nearness problem has a unique solution [4].
The matrix nearness problems have been extensively studied over the last thirty-two years. Much of the literature is available about some ad hoc mathematical methods that are not guaranteed to give the best solution to the problem. The method given by Knal and ten Burge [5] decomposes the matrix X as X = Y t Y and then minimizes the objective function over Y s columns. Lurie and Goldberg [6] minimize the quantity A R t R F 2 with R being as an upper triangular matrix having a unit 2-norm for it’s columns by using the Gauss–Newton method. In [7], Higham proposed an alternating projection algorithm using convex analysis, which has a linear convergence towards the solution of matrix nearness problems.
A more general matrix nearness problem than the above problem was studied in [8] by Malick, where he had used convex set and general linear constraints to replace positive semi-definiteness and unit d i a g ( X ) . Malick had applied the quasi-Newton method to dual problem after dualize the linear constraints on the d i a g ( X ) . A quadratically convergent Newton method for solving matrix nearness problem was studied by Qi and Sun [9], where they proposed a dual problem to matrix nearness problem and then proposed its solution. The numerical method developed by Toh [10] demands the solution of dense linear systems with dimension n 2 2 , mathematical construction of preconditioners and then apply them to compute the nearest correlation matrix X to solve the matrix nearness problem. A globally convergent Newton’s method proposed by Qi and Sun [9] computes the nearest correlation matrix X for the above matrix nearness problem.
In this work, we propose a low-rank ODE’s based method to compute a nearest correlation matrix X for a given n-dimensional matrix A. Our method works in two ways: first, it allows us to shift the smallest negative eigenvalue λ 1 of the perturbed matrix A + ϵ E , such that it becomes non-negative, that is, λ 1 0 . We compute the perturbation matrix E with d i a g ( E ) = 0 by constructing and then solving an associated optimization problem. Secondly, the construction and solution to optimization problem allow to shift all negative eigenvalues λ i , i = n 1 of given matrix A such that its spectrum becomes non-negative, that is, λ i 0 , i = 1 : n 1 .

Overview of the Article

In Section 2, we provide the preliminaries of our article. We give the definitions of symmetric, positive semi-definite matrix, a correlation matrix, and nearest correlation matrix.
Section 3 of our article is devoted to the problem under consideration. We briefly discuss the aim of matrix nearness problem. Section 4 of our article is dedicated to the computation of a gradient system of ordinary differential equations to localize the smallest negative eigenvalue λ 1 from the spectrum of an admissible perturbed matrix ( A + ϵ E ( t ) ) ) for given A R n , n .
The localization of eigenvalues λ 1 , λ 2 simultaneously from the spectrum of admissible perturbation matrix is addressed in Section 5. Finally, numerical experimentation and conclusion are presented in Section 6 and Section 7, respectively.

2. Preliminaries

Definition 1.
A matrix A R n , n is called symmetric if A t = A .
Definition 2.
A matrix A R n , n is called positive semi-definite if λ i ( A ) 0 , i = 1 : n , λ i ( A ) denotes the spectrum of matrix A .
Definition 3.
A square symmetric matrix A R n , n is called correlation matrix if ( i , j ) -entry is the correlation between the columns i and j of A.
Definition 4.
A matrix A R n , n is called a nearest correlation matrix if: A t = A , λ i ( A ) 0 , d i a g ( A ) = 1 , and a i j [ 1 , 1 ] , i j .

3. Matrix Nearness Problem

The matrix nearness problem aims the computation of a matrix X, which is symmetric, positive semi-definite, unit diagonal and having off diagonal entries between 1 and 1 for a given matrix A R n , n . The matrix A is not necessary to have it’s entries belonging to interval [ 1 , 1 ] . In fact, it can have real random entries. It is possible that all the eigenvalues of A are negative, or some of them are negative while others are positive. The negative eigenvalues are the basic hindrance to make A a nearest correlation matrix. Moreover, one must remain careful about the structure of A. The structure of A must have all properties, as defined in the preliminaries Section. For this purpose, we aim to compute the perturbation matrix E = E ( t ) , t having zero diagonal and a Frobenius norm bounded above by 1. In the very next Section, we present our methodology to compute the matrix E, which allows us to shift all negative eigenvalues of given A to make them non-negative.

4. Localization of Smallest Negative Eigenvalue λ 1

In this section, we aim to localize the smallest negative eigenvalue λ 1 of given matrix A R n , n . This involves shifting the smallest negative eigenvalue of perturbation matrix ( A + ϵ E ) with ϵ > 0 , a small positive parameter. The matrix E has structure with d i a g ( E ) = 0 . Moreover,
E ( t ) F = i j e i , j 2 1 , t .

4.1. Construction of Perturbation Matrix

We compute the perturbation matrix E ( t ) and then determine the direction Z = E ˙ ( t ) . The computation of E ˙ ( t ) indicates that how fast the smallest negative eigenvalue λ 1 ( t ) grows, that is, d d t ( λ i ( t ) ) > 0 . We need to compute the perturbed matrix ϵ E ( t ) with ϵ > 0 . For this purpose, we make use of eigenvalue problem
( A + ϵ E ( t ) ) η ( t ) = λ ( t ) η ( t ) ,
with η ( t ) being an eigenvector for smallest negative eigenvalue λ ( t ) . Furthermore, η ( t ) 2 1 .
In turn this implies
( η ( t ) * ( A + ϵ E ( t ) ) = λ ( t ) η ( t ) * .
By differentiating Equation (1) w.r.t. t, we have
( A + ϵ E ( t ) ) d d t η ( t ) + ϵ d d t ( E ( t ) ) η ( t ) = d d t ( λ ( t ) ) η ( t ) + λ ( t ) d d t ( η ( t ) ) .
Multiplying both sides with η ( t ) * ( t ) gives us
η ( t ) * ( A + ϵ E ( t ) ) d d t η ( t ) + ϵ η ( t ) * d d t ( E ( t ) ) η ( t ) = d d t ( λ ( t ) ) η ( t ) * η ( t ) + λ ( t ) η ( t ) * d d t ( η ( t ) ) .
As we know that,
η ( t ) * η ( t ) = η ( t ) , η ( t ) = η ( t ) 2 2 = 1 ,
This further implies that,
η ( t ) * ( A + ϵ E ( t ) ) d d t η ( t ) + ϵ η ( t ) * d d t ( E ( t ) ) η ( t ) = d d t ( λ ( t ) ) + λ ( t ) η ( t ) * + d d t ( η ( t ) ) .
Thus, we have
λ ( t ) η ( t ) * d d t ( η ( t ) ) + ϵ η ( t ) * d d t ( E ( t ) η ( t ) ) = d d t ( λ ( t ) ) + λ ( t ) η ( t ) * ) d d t ( η ( t ) ) .
Using Equation (2), we have that
λ ( t ) η ( t ) * = η ( t ) * ( A + ϵ E ( t ) ) .
Thus, Equation (3) becomes as
d d t ( λ ( t ) ) = ϵ η ( t ) * d d t ( E ( t ) ) η ( t ) .
We take η ( t ) * d d t ( η ( t ) ) = 0 in Equation (3), and Z = d d t ( E ( t ) ) = E ˙ ( t ) , results in the following optimization problem.

4.2. Formulation of Optimization Problem

The optimization problem given below allows one to determine the direction Z = E ˙ ( t ) so that the solution of the system of ODE’s indicates the maximum growth of the smallest negative eigenvalue λ ( t ) .
m a x ( η 1 * Z η 1 ) S u b j e c t t o Z , E ( t ) = 0 d i a g ( Z ) = 0 ,
with η 1 R n , 1 is eigenvector corresponding to smallest negative eigenvalue λ 1 . The notation ∗ represents complex conjugate transpose. The solution of the maximization problem addressed in Equation (5) is given as the following.

4.3. Lemma 4.2.1.

Consider that E ( t ) is a non-zero matrix with an admissible perturbation matrix
E ( t ) F i , j e i j 2 1 .
Let η 1 , η 1 * are non-zero eigen vectors associated with smallest negative eigenvalue λ 1 ( t ) . Then solution Z to maximization problem in Equation (5) is of the form
Z = P r o j ( η 1 η 1 * ) P r o j ( η 1 η 1 * , E ( t ) E ( t ) ,
having P r o j ( · ) as projection of direction matrix Z onto manifold of matrices E ( t ) .
Proof. 
The proof is based on the idea of computation of orthogonal projection on manifold of matrices, and for this we refer to [11]. □

4.4. A Gradient System of ODE’s

The solution matrix
Z = P r o j ( η 1 η 1 * ) P r o j ( η 1 η 1 * , E ( t ) E ( t ) ,
of the maximization problem addressed in Equation (5) allows following gradient system of ODEs on manifold of matrices E ( t ) ,
E ˙ ( t ) = P r o j ( η 1 η 1 * ) P r o j ( η 1 η 1 * , E ( t ) .

4.5. Characterization of Gradient System of ODE’s

The solution of gradient system of ODE’s addressed in Equation (7) possesses the following characteristic properties:
(i)
d d t ( λ 1 ( t ) ) > 0 ,
(ii)
E ˙ ( t ) = 0 d d t ( λ 1 ( t ) ) = 0 ,
(iii)
d d t ( λ 1 ( t ) ) = 0 E ( t ) P r o j ( η 1 η 1 * ) .

5. Localization of λ 1 ( t ) , λ 2 ( t )

Our goal in this Section is to transfer simultaneously the smallest negative eigenvalue λ 1 ( t ) and the negative eigenvalue λ 2 ( t ) which is next to λ 1 ( t ) from spectrum of perturbation matrix ( A + ϵ E ( t ) ) so that λ 1 ( t ) , λ 2 ( t ) becomes strictly positive.

5.1. Optimization Problem

The following maximization problem allows one to determine direction matrix Z = E ˙ ( t ) so that a solution associated with gradient system of ODEs obtained after solving the problem, which indicates maximum growth for λ 1 ( t ) and λ 2 ( t ) respectively. The maximization problem, in order to have a simultaneous maximum growth for both λ 1 ( t ) and λ 2 ( t ) , is given as
m a x ( η 1 * Z η 1 ) S u b j e c t t o η 2 * Z η 2 = η 1 * Z η 1 Z , E ( t ) = 0 d i a g ( Z ) = 0 .
Now, our aim is to give the solution corresponding to maximization problem addressed in Equation (8).

5.2. A Gradient System of ODE’s

The optimal solution to maximization addressed in Equation (8) is given by gradient system of ODE’s as
E ˙ ( t ) = ( 1 μ ) η 1 * η 1 + μ η 2 * η 2 μ { η 1 * η 1 η 2 * η 2 , E ( t ) η 1 * η 1 , E ( t ) } .
The gradient system of ODEs addressed in Equation (9) as a function K ( t ) with
K ( t ) = 1 a 11 ϵ 0 0 e 21 ( t ) 0 e n 1 ( t ) e n n ( t ) 1 a n n ϵ .
The choice of parameter ϵ could be sufficiently large enough so that d i a g ( A + ϵ E ( t ) ) = 1 and having
A = a 11 a 12 a 1 n a 21 a 22 a 2 n a n 1 a n 2 a n n ; E ( t ) = e 11 ( t ) e 12 ( t ) e 1 n ( t ) e 21 ( t ) e 22 ( t ) e 2 n ( t ) e n 1 ( t ) e n 2 ( t ) e n n ( t ) .
For parameter ϵ 0 ϵ , we have that ϵ 0 = d i a g ( A I ) ϵ . Next, we fix e 11 ( t ) , e 12 ( t ) , , e n n ( t ) and for parameter ϵ > > ϵ 0 , we obtain
A + ϵ E ( 0 ) = a 11 + ϵ e 11 ( t ) a 12 + ϵ e 12 ( t ) a 1 n + ϵ e 1 n ( t ) a 21 + ϵ e 21 ( t ) a 22 + ϵ e 22 ( t ) a 2 n + ϵ e 2 n ( t ) a n 1 + ϵ e n 1 ( t ) a n 2 + ϵ e n 2 ( t ) a n n + ϵ e n n ( t ) .
The admissible perturbation matrix A + ϵ E ( 0 ) ) takes the form
A + ϵ E ( 0 ) = 1 a 12 + ϵ e 12 ( t ) a 1 n + ϵ e 1 n ( t ) a 21 + ϵ e 21 ( t ) 1 a n n + ϵ e n n ( t ) a n 1 + ϵ e n 1 ( t ) a n n + ϵ e n n ( t ) 1
As, e 11 ( t ) = 1 a 11 ϵ , e 12 ( t ) = 1 a 22 ϵ , , e n n ( t ) = 1 a n n ϵ .
The initial perturbation matrix E ( 0 ) is given as
E ( 0 ) = 1 a 11 ϵ e 12 ( t ) e 1 n ( t ) e 21 ( t ) e n n ( t ) e n 1 ( t ) e n n ( t ) 1 a n n ϵ .
The initial perturbation matrix E ( 0 ) is decomposed into K T ( 0 ) and K ( 0 ) , the upper and lower triangular matrices respectively as,
K T ( 0 ) = 1 a 11 ϵ e 12 ( t ) e 1 n ( t ) 0 e n n ( t ) 0 0 1 a n n ϵ ; K ( 0 ) = 1 a 11 ϵ 0 0 e 21 ( t ) 0 e n 1 ( t ) e n n ( t ) 1 a n n ϵ .
Also,
2 K ( 0 ) 2 2 + i ( 1 a i i ) 2 ϵ 2 = 1 .
In a similar manner, the perturbation matrix E ( t ) has the structure
E ( t ) = 1 a 11 ϵ e 12 ( t ) e 1 n ( t ) e 21 ( t ) e n n ( t ) e n 1 ( t ) e n n ( t ) 1 a n n ϵ ,
and
K ( t ) = 1 2 1 i ( 1 a i i ) 2 ϵ 2 1 2 .
The maximization of eigenvalues λ 1 ( t ) and λ 2 ( t ) demands the computation of the perturbation matrix E ˙ ( t ) and projection of Z onto manifold of family of matrices E ( t ) . The change in the admissible perturbation matrix E ( t ) is given as
E ˙ ( t ) = 0 e ˙ 12 ( t ) e ˙ 1 n ( t ) e ˙ 21 ( t ) e ˙ n n ( t ) e ˙ n 1 ( t ) e ˙ n n ( t ) 0 .
Both K ˙ ( t ) and K ˙ T ( t ) are computed as
K ˙ ( t ) = 0 e ˙ 12 ( t ) e ˙ 1 n ( t ) 0 e ˙ n n ( t ) 0 0 0 ; K ˙ T ( t ) = 0 0 0 e ˙ 12 ( t ) 0 e ˙ 1 n ( t ) e ˙ n n ( t ) 0 .
The computation K ˙ ( t ) cause maximum growth of λ m i n ( t ) in the maximization problem as,
d d t ( λ m i n ( t ) ) = m a x ( η * E ˙ ( t ) η ) s u b j e c t t o E ( t ) , E ˙ ( t ) = 0 d i a g ( Z ) = 0 .
The solution to maximization problem addressed in Equation (12) is given by ODE E ˙ ( t ) as,
E ˙ ( t ) = P r o j ( η η * ) P r o j ( η η * , E ( t ) E ( t ) ,
having
P r o j ( η η * ) = η d i a g ( η 1 2 , η n 2 )
and has the form
P r o j ( η η * ) = 0 η 1 η 2 η 1 η n η 2 η 1 η n η n η n η 1 η n η n 0 .
The matrix P r o j ( η η * ) is decomposed into U ( t ) , an upper triangular matrix and lower triangular matrix L T ( t ) as
U ( t ) = 0 η 1 η 2 η 1 η n 0 η n η n 0 0 0 ,
L T ( t ) = 0 0 0 η 2 η 1 0 η n η 1 η n η n 0 .
Remark 1.
The upper triangular matrix U ( t ) = ( A + ϵ E ( t ) ) . The optimal solution to maximization problem addressed in Equation (8) has the form
E ˙ ( t ) = ( 1 μ ) P r o j ( η 1 η 1 * ) + μ P r o j ( η 2 η 2 * ) γ E ( t ) .
The solution E ˙ ( t ) addressed in Equation (13) is determined with the help of Euler’s method, that is,
E n + 1 = E n + h E ˙ n .
Finally, the eigenvalue problem
X = ( A + ϵ E ( t ) ) η ( t ) = λ ( t ) η ( t )
computes the non-negative spectrum for the perturbed system and it show that the matrix X is a nearest correlation matrix as it is symmetric, positive semi-definite, unit diagonal and off diagonal entries lies in the interval 1 and 1.

6. Numerical Experimentation

This section presents numerical experimentation for the computation of nearest correlation matrices for a given matrix A. For simplicity, we take A R n , n .
Example 1.
We consider a ten-dimensional matrix A as
A = 1 0 0 0 0 0 0 0 0 1 0 1 4 0 0 1 0 1 0 0 0 4 1 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 1 0 4 0 0 0 0 0 0 0 0 1 0 0 3 0 1 1 0 0 4 0 1 0 4 0 0 0 1 0 0 0 0 1 0 1 0 0 0 0 0 3 4 0 1 .
The eigenvalues of A are { 5.7249 , 2.8009 , 0.7529 , 0 , 1 , 1 , 2 , 2.9796 , 5.0768 , 7.2223 } which clearly contains the negative eigenvalues. The perturbation matrix E has a zero diagonal and help us to shift the negative eigenvalues of perturbed matrix ( A + ϵ E ) . The perturbation matrix E is computed as
E = 0 0.005 0.002 0.006 0.006 0.011 0.009 0.009 0.001 0.093 0.005 0 0.353 0.006 0.008 0.092 0.003 0.074 0.007 0.014 0.002 0.353 0 0.005 0.005 0.008 0.011 0.098 0.008 0.010 0.006 0.006 0.005 0 0.007 0.005 0.009 0.004 0.086 0.005 0.006 0.008 0.005 0.007 0 0.013 0.002 0.009 0.003 0.002 0.011 0.092 0.008 0.005 0.013 0 0.008 0.359 0.008 0.008 0.009 0.003 0.011 0.009 0.002 0.008 0 0.011 0.003 0.272 0.009 0.074 0.098 0.004 0.009 0.359 0.011 0 0.011 0.360 0.001 0.007 0.008 0.086 0.003 0.008 0.003 0.011 0 0.009 0.093 0.014 0.010 0.005 0.002 0.008 0.272 0.360 0.009 0 .
We compute matrix B 1 , which alters the eigenvalues of given matrix A; at the initial stage, this matrix also has three negative eigenvalues but different from those corresponding to A. The matrix B 1 and its eigenvalues are computed as
B 1 = 1 0.025 0.012 1.089 6.736 0.018 0.1993 0.154 0.826 0.792 0.025 1 2.901 0.017 0.072 0.918 0.048 0.893 0.014 0.080 0.012 2.901 1 0.013 0.362 0.194 0.068 0.963 0.001 0.026 1.089 0.017 0.013 1 2.888 0.095 0.312 0.004 0.017 0.074 6.736 0.072 0.362 2.888 1 0.665 1.759 0.982 2.213 0.104 0.018 0.918 0.194 0.095 0.665 1 0.360 3.120 0.034 0.153 0.199 0.048 0.068 0.312 1.759 0.360 1 0.096 0.205 2.209 0.154 0.893 0.963 0.004 0.982 3.120 0.096 1 0.011 3.360 0.826 0.014 0.001 0.017 2.213 0.034 0.205 0.011 1 0.039 0.792 0.080 0.026 0.074 0.104 0.153 2.209 3.360 0.039 1 .
The eigenvalues of B 1 are { 7.524 , 4.678 , 1.812 , 0.231 , 0.985 , 2.185 , 2.239 , 3.928 , 5.834 , 8.609 } and clearly is not a nearest correlation matrix. Matrix B 2 is also not a nearest correlation matrix but it has shifted the eigenvalues of A + ϵ E such that only one eigenvalue is negative from the spectrum. The matrix B 2 and its spectrum are given as
B 2 = 1 0.196 0.354 0.640 2.478 0.038 0.665 0.460 0.721 0.296 0.196 1 0.276 0.159 0.263 0.243 0.237 0.271 0.193 0.294 0.354 0.276 1 0.218 0.435 0.243 0.312 0.326 0.277 0.286 0.640 0.159 0.218 1 1.215 0.052 0.574 0.369 0.464 0.203 2.478 0.263 0.435 1.215 1 0.001 0.961 0.632 1.092 0.308 0.038 0.243 0.243 0.052 0.001 1 0.007 0.230 0.050 0.267 0.665 0.237 0.312 0.574 0.961 0.007 1 0.285 0.471 0.106 0.460 0.271 0.326 0.369 0.632 0.230 0.285 1 0.344 0.264 0.721 0.193 0.277 0.464 1.092 0.050 0.471 0.344 1 0.128 0.296 0.294 0.286 0.203 0.308 0.267 0.106 0.264 0.128 1 .
{ 3.560 , 0.289 , 1.195 , 1.251 , 1.307 , 1.335 , 1.395 , 1.570 , 1.626 } respectively. The matrix B 3 is a nearest correlation matrix, as not only does it have a unit diagonal, but its symmetric, positive semi-definite and off diagonal entries lie within the interval [ 1 , 1 ] . The matrix B 3 and its positive spectrum are given as follows.
B 3 = 1 0.054 0.111 0.015 0.036 0.097 0.057 0.113 0.069 0.104 0.054 1 0.174 0.090 0.094 0.125 0.109 0.094 0.065 0.115 0.111 0.174 1 0.124 0.105 0.073 0.102 0.140 0.033 0.134 0.015 0.090 0.124 1 0.050 0.061 0.040 0.078 0.079 0.049 0.036 0.094 0.105 0.050 1 0.055 0.099 0.089 0.034 0.112 0.097 0.125 0.073 0.061 0.055 1 0.156 0.236 0.089 0.118 0.057 0.109 0.102 0.040 0.099 0.156 1 0.130 0.105 0.158 0.113 0.094 0.140 0.078 0.089 0.236 0.130 1 0.096 0.266 0.069 0.065 0.033 0.079 0.034 0.089 0.105 0.096 1 0.094 0.104 0.115 0.134 0.049 0.112 0.118 0.158 0.266 0.094 1 .
The positive spectrum of nearest correlation matrix is {0.061, 0.936, 1.004, 1.025, 1.067, 1.084, 1.131, 1.153, 1.217, 1.319}.
Example 2.
We consider an eight dimensional matrix A as
A = 2 1 3 0 2 1 0 2 1 2 4 1 3 1 2 1 3 4 2 1 2 0 2 1 0 1 1 2 1 0 4 2 2 3 2 1 2 1 3 0 1 1 0 0 1 2 1 4 0 2 2 4 3 1 2 1 2 1 1 2 0 4 1 2 .
The perturbation matrix E has a zero diagonal and helps us to shift the negative eigenvalues of perturbed matrix ( A + ϵ E ) . The perturbation matrix E is computed as
E = 0 0.073 0.199 0.003 0.136 0.066 0.001 0.131 0.073 0 0.273 0.070 0.206 0.066 0.130 0.065 0.199 0.273 0 0.073 0.131 0.001 0.129 0.073 0.003 0.070 0.073 0 0.070 0.001 0.261 0.130 0.136 0.206 0.131 0.070 0 0.071 0.206 0 0.066 0.066 0.001 0.001 0.071 0 0.063 0.275 0.001 0.130 0.129 0.261 0.206 0.063 0 0.063 0.131 0.065 0.073 0.130 0 0.275 0.063 0 .
We compute matrix B 1 , which alters the eigenvalues of given matrix A; at the initial stage, this matrix also has three negative eigenvalues but different from those corresponding to A. The matrix B 1 and its eigenvalues are computed as
B 1 = 1 0.846 3.158 0.876 0.326 0.875 0.045 0.8764 0.846 1 0.418 1.143 2.025 1.677 1.050 0.792 3.158 0.418 1 0.543 1.252 0.826 2.071 0.678 0.876 1.143 0.543 1 0.631 0.718 2.403 1.686 0.326 2.025 1.252 0.631 1 1.631 0.372 1.696 0.875 1.677 0.826 0.718 1.631 1 1.783 2.473 0.045 1.050 2.071 2.403 0.372 1.783 1 2.224 0.876 0.792 0.678 1.686 1.696 2.473 2.224 1 .
The eigenvalues of B 1 are { 5.691 , 3.971 , 0.152 , 0.540 , 3.484 , 3.578 , 4.002 , 6.208 } and clearly is not a nearest correlation matrix. Matrix B 2 is also not a nearest correlation matrix but it has shifted the eigenvalues of A + ϵ E , such that two eigenvalues are negative from the spectrum. The matrix B 2 and its spectrum are given as
B 2 = 1 1.114 1.683 0.370 0.020 0.717 0.208 0.849 1.114 1 0.363 0.807 2.048 1.146 0.921 0.668 1.683 0.363 1 0.451 1.263 0.679 1.089 0.348 0.370 0.807 0.451 1 0.885 1.044 1.469 1.275 0.020 2.048 1.263 0.885 1 1.031 0.810 1.2782 0.717 1.146 0.679 1.044 1.031 1 1.198 1.592 0.208 0.921 1.089 1.469 0.810 1.198 1 1.573 0.849 0.668 0.348 1.275 1.278 1.592 1.573 1 .
{ 4.211 , 2.776 , 0.745 , 1.701 , 2.316 , 2.514 , 2.879 , 4.832 } respectively. The matrix B 3 is a nearest correlation matrix, as not only does it have a unit diagonal, but its symmetric, positive semi-definite and off diagonal entries lie within the interval [ 1 , 1 ] . The matrix B 3 and its positive spectrum are given as follows.
B 3 = 1 0.013 0.142 0.006 0.047 0.067 0.031 0.139 0.013 1 0.087 0.103 0.113 0.055 0.011 0.150 0.142 0.087 1 0.099 0.018 0.042 0.116 0.014 0.006 0.103 0.099 1 0.011 0.042 0.116 0.032 0.047 0.113 0.018 0.011 1 0.055 0.152 0.025 0.067 0.055 0.042 0.042 0.055 1 0.063 0.071 0.031 0.011 0.116 0.116 0.152 0.063 1 0.025 0.139 0.150 0.014 0.032 0.025 0.071 0.025 1 .
The spectrum of nearest correlation matrix is obtained as { 0.680 , 0.790 , 0.855 , 0.965 , 1.045 , 1.081 , 1.188 , 1.393 } .

7. Conclusions

In this article, we have presented a low-rank ordinary differential equations-based technique to compute the nearest correlation matrix; that is, symmetric, positive semidefinite, unit diagonal and off-diagonal entries between 1 and 1 for a given n-dimensional real-valued matrix. The proposed methodology is based on transforming the negative spectrum of the original matrix, which involves the computation of perturbation level and an admissible perturbation having zero diagonal. The numerical experimentation turns over the desired perturbations for the computation of nearest correlation matrix for randomly chosen real-valued matrices.

Author Contributions

M.-U.R. contributed in conceptualization, methodology, validation of results, while J.A. and K.A. contributed in writing-review and editing. All authors have read and agreed to the published version of the manuscript.

Funding

J. Alzabut and K. Abodayeh would like to thank Prince Sultan University for funding this work through research group Nonlinear Analysis Methods in Applied Mathematics (NAMAM) group number RG-DES-2017-01-17.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Finger, C. A methodology to stress correlations. RiskMetrics Monit. 1997, 4, 3–11. [Google Scholar]
  2. Higham, N.J. Computing a nearest symmetric positive semidefinite matrix. Linear Algebra Its Appl. 1988, 103, 103–118. [Google Scholar] [CrossRef]
  3. Saygun, T.; Epperlein, E.; Christofides, N. Correlation stress testing for value-at-risk. J. Risk 2003, 5, 75–89. [Google Scholar]
  4. Deutsch, F. Best Approximation. In Best Approximation in Inner Product Spaces; Springer: New York, NY, USA, 2001; pp. 21–32. [Google Scholar]
  5. Knol, D.L.; ten Berge, J. Least-squares approximation of an improper correlation matrix by a proper one. Psychometrika 1989, 54, 53–61. [Google Scholar] [CrossRef]
  6. Lurie, P.M.; Goldberg, M.S. An approximate method for sampling correlated random variables from partially-specified distributions. Manag. Sci. 1998, 44, 203–218. [Google Scholar] [CrossRef]
  7. Higha, N.J. Computing the nearest correlation matrix—A problem from finance. IMA J. Numer. Anal. 2002, 22, 329–343. [Google Scholar] [CrossRef] [Green Version]
  8. Malick, J. A dual approach to semidefinite least-squares problems. SIAM J. Matrix Anal. Appl. 2004, 26, 272–284. [Google Scholar] [CrossRef] [Green Version]
  9. Qi, H.; Sun, D. A quadratically convergent Newton method for computing the nearest correlation matrix. SIAM J. Matrix Anal. Appl. 2006, 28, 360–385. [Google Scholar] [CrossRef]
  10. Toh, K.-C. An inexact primal–dual path following algorithm for convex quadratic SDP. Math. Program. 2008, 112, 221–254. [Google Scholar] [CrossRef]
  11. Guglielmi, N.; Manetta, M. Approximating real stability radii. IMA J. Numer. Anal. 2015, 35, 1402–1425. [Google Scholar] [CrossRef]
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Rehman, M.-U.; Alzabut, J.; Abodayeh, K. Computing Nearest Correlation Matrix via Low-Rank ODE’s Based Technique. Symmetry 2020, 12, 1824. https://doi.org/10.3390/sym12111824

AMA Style

Rehman M-U, Alzabut J, Abodayeh K. Computing Nearest Correlation Matrix via Low-Rank ODE’s Based Technique. Symmetry. 2020; 12(11):1824. https://doi.org/10.3390/sym12111824

Chicago/Turabian Style

Rehman, Mutti-Ur, Jehad Alzabut, and Kamaleldin Abodayeh. 2020. "Computing Nearest Correlation Matrix via Low-Rank ODE’s Based Technique" Symmetry 12, no. 11: 1824. https://doi.org/10.3390/sym12111824

APA Style

Rehman, M. -U., Alzabut, J., & Abodayeh, K. (2020). Computing Nearest Correlation Matrix via Low-Rank ODE’s Based Technique. Symmetry, 12(11), 1824. https://doi.org/10.3390/sym12111824

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop