# Error with Bayesian Gaussian Process "PosDefException: matrix is not positive definite; Cholesky factorization failed."

**URL:** https://discourse.julialang.org/t/error-with-bayesian-gaussian-process-posdefexception-matrix-is-not-positive-definite-cholesky-factorization-failed/85683
**Category:** Statistics
**Created:** [August 12, 2022, 9:15pm UTC](https://discourse.julialang.org/t/error-with-bayesian-gaussian-process-posdefexception-matrix-is-not-positive-definite-cholesky-factorization-failed/85683 "2022-08-12T21:15:13Z")
**Posts on this page:** 3
**Page:** 1

<div class="post-metadata">

### Author: ![vickiluv](https://avatars.discourse-cdn.com/v4/letter/v/a183cd/32.png) [@vickiluv](https://discourse.julialang.org/u/vickiluv)
#### Post date: [August 12, 2022, 9:15pm UTC](https://discourse.julialang.org/t/error-with-bayesian-gaussian-process-posdefexception-matrix-is-not-positive-definite-cholesky-factorization-failed/85683/1 "2022-08-12T21:15:13Z")

</div>

I got the error of “PosDefException” when trying Bayesian Gaussian Process in Turing.jl. The Bayesian framework looks as follows

```julia
@model function GP_μ(μ, D; cov_fn::Function=exp_cov_fn)
    n_station = length(μ)

    # weakly informative prior over mean field for the μ_i
    μ_mean ~ Normal(100, 25)
	
    # weakly informative priors over the kernel variance parameter
    # kernel variance parameter is > 0
    μ_kernelvariance ~ InverseGamma(10, 50) # 10 ± 5 roughly

    # weakly informative priors over kernel length scale ℓ
    min_dist = minimum(D) + 1E-3
    max_dist = maximum(D)
    μ_kernellength ~ Uniform(min_dist, max_dist)
	
    ########################
    # GAUSSIAN PROCESS MODEL
    ########################

    # the kernels are computed exactly
    μ_kernel = cov_fn(D, μ_kernelvariance_log, μ_kernellength)
	
    # process layer implements the Gaussian Process
    δ = 1E-6 # add an offset for numerical stability
    μ ~ MvNormal(μ_mean * ones(n_station), Matrix(Hermitian(μ_kernel + δ * I)))

end

function calc_μ_true(lon, lat)
    dist = Euclidean()(CENTER, [lon, lat])
    return 100 + 25 * sin(2 * π * dist / 5 + 3π / 2) + rand(Normal(0, 10))/5
end

n_station = 30 # number of synthetic stations
n_obs = 30 # number of obs per station

# get the true values
lons = rand(Uniform(LON_RANGE[1], LON_RANGE[2]), n_station)
lats = rand(Uniform(LAT_RANGE[1], LAT_RANGE[2]), n_station)
μ_true = calc_μ_true.(lons, lats)
distance = calc_dist(hcat(lats, lons))

model_μ = GP_μ(μ_true, distance)
chains = begin
	mcmc_sampler = DynamicNUTS()
	n_per_chain = 1000
	nchains = 2
	sample(model_μ, mcmc_sampler, MCMCThreads(), n_per_chain, nchains; drop_warmup=true)
end

```

I tried adding a noise term, using Hermitian(), Symmetric(), but still got the error of either “matrix is not positive definite” or “matrix is not Hermitian” for the covariance I suppose. Any idea how I might solve this?

---

<div class="post-metadata">

### Author: ![blackeneth](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/blackeneth/32/10353_2.png) [@blackeneth](https://discourse.julialang.org/u/blackeneth)
#### Post date: [August 12, 2022, 10:56pm UTC](https://discourse.julialang.org/t/error-with-bayesian-gaussian-process-posdefexception-matrix-is-not-positive-definite-cholesky-factorization-failed/85683/2 "2022-08-12T22:56:44Z")

</div>

What do you get when you just compute the cholesky factorization using cholesky()?

---

<div class="post-metadata">

### Author: ![vickiluv](https://avatars.discourse-cdn.com/v4/letter/v/a183cd/32.png) [@vickiluv](https://discourse.julialang.org/u/vickiluv)
#### Post date: [August 13, 2022, 3:57am UTC](https://discourse.julialang.org/t/error-with-bayesian-gaussian-process-posdefexception-matrix-is-not-positive-definite-cholesky-factorization-failed/85683/4 "2022-08-13T03:57:04Z")

</div>

I got results like this

```julia
LinearAlgebra.Cholesky{Float64, Matrix{Float64}}([7.712549336531022 1.847900191359362 … 1.4450590383187025 0.251966480687946; 14.252021394844194 7.4879023865966055 … 0.0032314205728996553 0.5100249881619099; … ; 11.145089127233064 2.694521435254619 … 4.61337750842342 -1.0599196726807111e-5; 1.9433039134578745 4.2846262339608705 … 0.45873498713978184 4.399661715685264], 'U', 0)

```

I got the PosDefException error from bayesian sampling, but computing cholesky() or isposdef() on the MvNormal covariance (`μ_kernel`) seemed fine.