# Verify a matrix is positive semi-definite

**URL:** https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329
**Category:** General Usage
**Created:** [April 20, 2019, 2:28am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329 "2019-04-20T02:28:49Z")
**Posts on this page:** 14
**Page:** 1

<div class="post-metadata">

### Author: ![Mr.Robot](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mr.robot/32/8052_2.png) [@Mr.Robot](https://discourse.julialang.org/u/Mr.Robot)
#### Post date: [April 20, 2019, 2:28am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/1 "2019-04-20T02:28:49Z")

</div>

I am trying to numerically verify that

> A symmetric matrix \mathbf{A} is positive semidefinite if and only if it is a covariance matrix.

Then I need to verify in both directions, i.e.

1. Given a positive semidefinite matrix \mathbf{A}, show that it is a covariance matrix.
2. Given a covariance matrix, show that it is positive semidefinite.

However, I am not sure

1. What properties should a matrix have to be a covariance matrix.
2. I know I could generate a covariance matrix using the following and I know that `cov` is positive semidefinite if and only if all of its eigenvalues are non-negative. But it turns out that `minimum(eigvals(cov))` is a negative number close to 0 (on the order \sim 10^{-15}), I am not sure if I could conclude that `cov` is positive semidefinite since numerical reasons.

```julia
n = 100
u = randn(n);
cov = u * u'

```

Any input will be appreciated.

---

<div class="post-metadata">

### Author: ![aharoun](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/aharoun/32/6887_2.png) [@aharoun](https://discourse.julialang.org/u/aharoun)
#### Post date: [April 20, 2019, 5:14am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/3 "2019-04-20T05:14:22Z")

</div>

The issue is that what you are creating is not a covariance matrix. Consider this

```julia
julia> n = 100

julia> x = randn(n,3); # 3 random variables with sample size n

julia> demean_x = x .- mean(x,dims=1)

julia> cov_x = demean_x'*demean_x./(n-1) # or simply use cov function from Statistics

julia> minimum(eigvals(cov_x))>=0
true

```

Julia also has a function `isposdef` under `LinearAlgebra`, it checks if the matrix is **positive definite**.

```julia
julia> using LinearAlgebra

julia> isposdef(cov_x)
true

```

---

<div class="post-metadata">

### Author: ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)
#### Post date: [April 20, 2019, 5:23am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/4 "2019-04-20T05:23:25Z")

</div>

> [@Mr.Robot](#):
>
> Then I need to verify in both directions, i.e.
> 
> 1. Given a positive semidefinite matrix A\mathbf{A}, show that it is a covariance matrix.
> 2. Given a covariance matrix, show that it is positive semidefinite.

I am not sure I understand the whole context, but this is a [well-known property](https://en.wikipedia.org/wiki/Definiteness_of_a_matrix#Connections) (covariance matrices are psd, and psd matrices are covariance matrices). I don’t know how you could verify a statement about an uncountably infinite set numerically.

---

<div class="post-metadata">

### Author: ![Mr.Robot](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mr.robot/32/8052_2.png) [@Mr.Robot](https://discourse.julialang.org/u/Mr.Robot)
#### Post date: [April 20, 2019, 6:04am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/5 "2019-04-20T06:04:27Z")

</div>

It is actually a homework question. I think the point of this question is that professor wants us to get familiar with basic Julia programming.

---

<div class="post-metadata">

### Author: ![mcreel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcreel/32/30088_2.png) [@mcreel](https://discourse.julialang.org/u/mcreel)
#### Post date: [April 20, 2019, 6:43am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/6 "2019-04-20T06:43:09Z")

</div>

To verify numerically 1, you may find [https://github.com/mcreel/Econometrics/blob/master/Examples/GLS/cholesky.jl](https://github.com/mcreel/Econometrics/blob/master/Examples/GLS/cholesky.jl) useful.

To verify 2, just check the eigenvalues. Your code is missing the step of averaging to compute the sample covariance, and is not keeping the dimensionality fixed.

---

<div class="post-metadata">

### Author: ![tkoolen](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkoolen/32/1603_2.png) [@tkoolen](https://discourse.julialang.org/u/tkoolen)
#### Post date: [April 20, 2019, 4:31pm UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/7 "2019-04-20T16:31:47Z")

</div>

A Cholesky decomposition requires positive definiteness though; semidefiniteness is not enough.

---

<div class="post-metadata">

### Author: ![Rouzaire](https://avatars.discourse-cdn.com/v4/letter/r/85f322/32.png) [@Rouzaire](https://discourse.julialang.org/u/Rouzaire)
#### Post date: [April 9, 2020, 8:23am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/8 "2020-04-09T08:23:26Z")

</div>

Sorry for reviving this topic but I think that my question is quite similar to the title of this one. Let me know if I should have created a new one.

I have a Gram matrix **A** constructed from a Matérn covariance, so mathematically **A** is posdef. The thing is that for small such matrices, `posdef`( **A** ) is true, but for large such matrices, `posdef`( **A** ) becomes false. I guess that this arises from rounding errors but do you know any way to deal with this numeric issue ?

---

<div class="post-metadata">

### Author: ![RoyiAvital](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/royiavital/32/571_2.png) [@RoyiAvital](https://discourse.julialang.org/u/RoyiAvital)
#### Post date: [April 9, 2020, 8:30am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/9 "2020-04-09T08:30:04Z")

</div>

This is probably due to floating point calculations.  
Just add to it scaled identity matrix. Scaling factor should have the value of the minimum eigen value of your matrix (Maybe additional `eps()`).

---

<div class="post-metadata">

### Author: ![Rouzaire](https://avatars.discourse-cdn.com/v4/letter/r/85f322/32.png) [@Rouzaire](https://discourse.julialang.org/u/Rouzaire)
#### Post date: [April 9, 2020, 8:35am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/10 "2020-04-09T08:35:32Z")

</div>

Oh, that’s the so-called _jitter_, right ? Let’s say that one obtains a matrix **Ã** = **A** + tau^2\ ***I**  
But then, does a random number generation `rand(MvNormal(` **0** , **Ã** `))` have the same properties than `rand(MvNormal(` **0** , **A** `))` ?

---

<div class="post-metadata">

### Author: ![RoyiAvital](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/royiavital/32/571_2.png) [@RoyiAvital](https://discourse.julialang.org/u/RoyiAvital)
#### Post date: [April 9, 2020, 8:39am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/11 "2020-04-09T08:39:17Z")

</div>

Well, since A isn’t a PSD it is not a valid covariance matrix, so the questions is not well defined :-).  
But generally speaking, yes, it will be quite similar.

Another method is to threshold negative eigen values.  
So you calculate the Eigen Decomposition of A , and any negative value you set to 0 and then rebuild the matrix. It might be closer to what you need.

---

<div class="post-metadata">

### Author: ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)
#### Post date: [April 9, 2020, 9:33am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/12 "2020-04-09T09:33:37Z")

</div>

> [@RoyiAvital](#):
>
> set to 0 and then rebuild the matrix

Rebuilding may not be necessary, if factors are needed in the later steps anyway. Which is the case for multivariate normals.

Generally, if the problem permits, it is best not to form A at all, but work with a representation that either precludes impossible values, or makes it very easy to regularize them. I am not familiar with Matérn covariance, but for plain vanilla covariance one could just proceed from the QR decomposition of the data.

---

<div class="post-metadata">

### Author: ![Rouzaire](https://avatars.discourse-cdn.com/v4/letter/r/85f322/32.png) [@Rouzaire](https://discourse.julialang.org/u/Rouzaire)
#### Post date: [April 13, 2020, 8:56am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/13 "2020-04-13T08:56:44Z")

</div>

I come back to say that was a perfect solution to solve my problem, thank you !

---

<div class="post-metadata">

### Author: ![RoyiAvital](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/royiavital/32/571_2.png) [@RoyiAvital](https://discourse.julialang.org/u/RoyiAvital)
#### Post date: [April 13, 2020, 9:07am UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/14 "2020-04-13T09:07:21Z")

</div>

I am happy to hear I could assist.

Enjoy…

---

<div class="post-metadata">

### Author: ![baggepinnen](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/baggepinnen/32/693_2.png) [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)
#### Post date: [April 13, 2020, 12:04pm UTC](https://discourse.julialang.org/t/verify-a-matrix-is-positive-semi-definite/23329/15 "2020-04-13T12:04:00Z")

</div>

See also [GitHub - timholy/PositiveFactorizations.jl: Positive-definite "approximations" to matrices](https://github.com/timholy/PositiveFactorizations.jl)
