# Degenerate Multivariate Normal Sampling

**URL:** https://discourse.julialang.org/t/degenerate-multivariate-normal-sampling/103825
**Category:** Statistics
**Created:** [September 13, 2023, 7:09pm UTC](https://discourse.julialang.org/t/degenerate-multivariate-normal-sampling/103825 "2023-09-13T19:09:40Z")
**Posts on this page:** 2
**Page:** 1

<div class="post-metadata">

### Author: ![Jad\_Zeitouni](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jad_zeitouni/32/50678_2.png) [@Jad\_Zeitouni](https://discourse.julialang.org/u/Jad_Zeitouni)
#### Post date: [September 13, 2023, 7:09pm UTC](https://discourse.julialang.org/t/degenerate-multivariate-normal-sampling/103825/1 "2023-09-13T19:09:40Z")

</div>

It seems this topic keeps coming up over and over again (see [here](https://discourse.julialang.org/t/multivariate-normal-with-positive-semi-definite-covariance-matrix/3029) and [here](https://github.com/JuliaStats/Distributions.jl/issues/1219)), with a lot of [efforts to resolve](https://github.com/JuliaStats/Distributions.jl/pull/1552/f) it but it gets stuck.

My code is reliant on using `MvNormal` for sampling. I could write my own code that checks if the covariance matrix I am using is positive definite or positive semi-definite. Then create my own special constructor for `MvNormal` in the degenerate case and my own `rand` method, but I feel very uneasy about overwriting a package.

Is there a cleaner way to do this?

---

<div class="post-metadata">

### Author: ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)
#### Post date: [October 28, 2023, 2:22pm UTC](https://discourse.julialang.org/t/degenerate-multivariate-normal-sampling/103825/2 "2023-10-28T14:22:01Z")

</div>

Basically it means that your multivariate gaussian random vector is a linear transformation of a lower-dimensional standard gaussian random vector.

To find the right linear transformation, remove first the mean and then compute a (r,d) square root Matrix of your covariance. There are several algorithms to do that, and if the one used by MvNormal is not working in this particular case then pick another one by googling “square root of semi-positive deginite matrix in Julia”.

Then use `X = mu + Gamma * MvNormal(r)`, where mu is the mean, gamma the matrix square root, and r the rank.
