# Cholesky\_factor\_corr initializing with invalid values

**URL:** https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167
**Category:** Developers
**Created:** [February 1, 2018, 12:58am UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167 "2018-02-01T00:58:25Z")
**Posts on this page:** 9
**Page:** 1

<div class="post-metadata">

### Author: ![JustinYap](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/justinyap/32/664_2.png) [@JustinYap](https://discourse.mc-stan.org/u/JustinYap)
#### Post date: [February 1, 2018, 12:58am UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/1 "2018-02-01T00:58:25Z")

</div>

Thanks to those who have answered my previous questions, it is great that there is an active and helpful community for Stan.

Today I have potentially what might be a bug with `cholesky_factor_corr`. Basically when the size of the correlation matrix gets large, e.g. 40 x 40, I run into issues with the initial values for Cholesky factor correlation parameters:

```
Rejecting initial value:
  Error evaluating the log probability at the initial value.
Exception: multi_normal_lpdf: LDLT_Factor of covariance parameter is not positive definite. last conditional variance is 6.66134e-015. (in 'model3d5066d27a9a_0b961540b66c542856a20ec83237a5e1' at line 18)

```

I have a simple example in R below:

```
stan.code <- "
data {
    int K;
    int R;
    vector[K] draws[R];
    vector[K] sigma;
    vector[K] mu_0;
}
parameters {
    cholesky_factor_corr[K] L_omega;
}
transformed parameters {
    matrix[K, K] Sigma;
    Sigma = quad_form_diag(tcrossprod(L_omega), sigma);
}
model {
    draws ~ multi_normal(mu_0, Sigma);
}
"

n.draws <- 200
n.variables <- 40

draws <- matrix(rnorm(n.variables * n.draws),
                nrow = n.draws,
                ncol = n.variables)

stan.data <- list(
    K = n.variables,
    R = n.draws,
    draws = draws,
    sigma = rep(1, n.variables),
    mu_0 = rep(0, n.variables)
)

fit <- stan(model_code = stan.code,
            data = stan.data,
            iter = 100,
            chains = 1)

```

If `n.variables` is reduced to 20 or I set `init = 0` then everything works fine, so it appears to be an issue with initialization of the Cholesky factor correlation parameter when the matrix is large. Right now my solution is to pass in an identity matrix as the initial value. Has anyone reported this before?

---

<div class="post-metadata">

### Author: ![bgoodri](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bgoodri/32/4451_2.png) [@bgoodri](https://discourse.mc-stan.org/u/bgoodri)
#### Post date: [February 1, 2018, 1:24am UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/2 "2018-02-01T01:24:04Z")

</div>

> [@JustinYap](#):
>
> Has anyone reported this before?

Yes, but it is not really a bug, just a numerical problem. You could also do `init_r = ` some number less than the default of 2.

---

<div class="post-metadata">

### Author: ![JustinYap](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/justinyap/32/664_2.png) [@JustinYap](https://discourse.mc-stan.org/u/JustinYap)
#### Post date: [February 1, 2018, 2:00am UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/3 "2018-02-01T02:00:41Z")

</div>

Thanks, but that would affect all parameters and the number to use needs to be found with trial and error and may not work with different data. Are there plans to fix this numerical problem? Also, what are downsides of setting `init = 0` all the time, apart from not being able to randomize the starting point?

---

<div class="post-metadata">

### Author: ![bgoodri](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bgoodri/32/4451_2.png) [@bgoodri](https://discourse.mc-stan.org/u/bgoodri)
#### Post date: [February 1, 2018, 2:26am UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/4 "2018-02-01T02:26:48Z")

</div>

> [@JustinYap](#):
>
> but that would affect all parameters and the number to use needs to be found with trial and error and may not work with different data

The unconstrained parameters are all initialized on [-2,2], which if anything is too big.

> [@JustinYap](#):
>
> Are there plans to fix this numerical problem?

No, as it is not really a problem with the code but a fundamental limitation of how computers do math.

> [@JustinYap](#):
>
> what are downsides of setting init = 0 all the time, apart from not being able to randomize the starting point?

For this model, that is about it. If you do anything with `unit_vector`s, the sampling won’t start if `init = 0`.

---

<div class="post-metadata">

### Author: ![Bob\_Carpenter](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bob_carpenter/32/9230_2.png) [@Bob\_Carpenter](https://discourse.mc-stan.org/u/Bob_Carpenter)
#### Post date: [February 7, 2018, 10:56pm UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/5 "2018-02-07T22:56:54Z")

</div>

> [@JustinYap](#):
>
> Are there plans to fix this numerical problem?

> [@bgoodri](#):
>
> No, as it is not really a problem with the code but a fundamental limitation of how computers do math.

Sometimes these problems can be fixed.

What’s failing is the Eigen library impelemntation of Cholesky factorization for badly conditioned matrices. The bad conditioning comes up from our initializations. We don’t know how to initialize robustly in general as it depends on the scale of the problem, which we can’t really know in advance.

Our code’s not well set up to change initialization for a single variable type just for initialization. We could tinker with trying to make the transforms more robust, but we don’t really know how to do that—it’s an open problem.

---

<div class="post-metadata">

### Author: ![bgoodri](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bgoodri/32/4451_2.png) [@bgoodri](https://discourse.mc-stan.org/u/bgoodri)
#### Post date: [February 8, 2018, 12:07am UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/6 "2018-02-08T00:07:10Z")

</div>

Or the diagonal of the Cholesky factor underflows to zero and Eigen is doing the right thing by saying it is not positive definite.

---

<div class="post-metadata">

### Author: ![Bob\_Carpenter](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bob_carpenter/32/9230_2.png) [@Bob\_Carpenter](https://discourse.mc-stan.org/u/Bob_Carpenter)
#### Post date: [February 8, 2018, 12:18am UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/7 "2018-02-08T00:18:13Z")

</div>

Right. If the true value can’t be represented in floating point without underflowing, Eigen has no choice.

If it happens on an intermediate calculation, there’s then the issue of trying to fix it at the possible cost of speed for non-boundary cases and increased doc and maintenance burden.

---

<div class="post-metadata">

### Author: ![kn2465](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/kn2465/32/11528_2.png) [@kn2465](https://discourse.mc-stan.org/u/kn2465)
#### Post date: [July 13, 2021, 6:20pm UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/8 "2021-07-13T18:20:49Z")

</div>

Hi there, I’m facing the same issue when initializing a Kronecker Gaussian Process model, wonder if there has been any update to this open problem yet?

---

<div class="post-metadata">

### Author: ![andrjohns](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/andrjohns/32/15297_2.png) [@andrjohns](https://discourse.mc-stan.org/u/andrjohns)
#### Post date: [July 14, 2021, 12:19am UTC](https://discourse.mc-stan.org/t/cholesky-factor-corr-initializing-with-invalid-values/3167/9 "2021-07-14T00:19:36Z")

</div>

One approach is to try initialising using an identity matrix, that’s worked well for me in the past
