ENH: stats.gaussian_kde: replace use of inv_cov in logpdf

[mdhaber]

Reference issue

gh-16692

What does this implement/fix?

gh-16692 used Cholesky decomposition to avoid the inversion of the covariance matrix in gaussian_kde.pdf.
We wanted to wait for gh-15493 to finish before implementing that change in gaussian_kde.logpdf.
Now that gh-15493 is done, this makes the change to gaussian_kde.logpdf.

It also simplifies what we did in gh-16692. In that, we found a way to use Cholesky decomposition to compute $L^T x$, where $L L^T = \Sigma^{-1}$ (that is, $L$ is a Cholesky factor of the precision matrix / inverse covariance matrix).
Ultimately, we don't need $L^T x$; it was just one way to get to $x^T \Sigma^{-1} x$. There is a simpler way to get that while still avoiding the matrix inversion: $x^T \Sigma^{-1} x = y^T y$, where $C y = x$ and $CC^T = \Sigma$ (that is, $C$ is a Cholesky factor of the original covariance matrix rather than the precision matrix). The short calculation is shown here.

@steppi I think someone else can review this since it's simpler than gh-16692, but I thought you might find it interesting that those permutations weren't necessary to get the end result.

[mdhaber]

Failures appear to be unrelated.

[steppi]

@steppi I think someone else can review this since it's simpler than gh-16692, but I thought you might find it interesting that those permutations weren't necessary to get the end result.

Oh right, that makes perfect sense. If you start with the cholesky decomposition of the covariance matrix then you end up with something that isn’t quite a cholesky decomposition of the precision matrix because the left factor is upper triangular instead of lower. This factorization is still workable though because you can just call solve_triangular with lower=True in the end.

[steppi]

Looks good to me.

[mdhaber]

Thanks @steppi.

@jjerphan does this look good to you from a Cython perspective? To avoid too many changes at once, I didn't release the GIL. That can be a separate PR, if desired.

[jjerphan]

@jjerphan does this look good to you from a Cython perspective? To avoid too many changes at once, I didn't release the GIL. That can be a separate PR, if desired.

This LGTM from a Cython perspective -- I haven't looked at the maths, yet. Doing another PR for improving this implementation w.r.t Cython technicalities is sensible to me. 👍

Let me know if you need another review for this PR generally.

[mdhaber]

@steppi would you like @jjerphan to review the math, too, or are you comfortable merging?

[steppi]

@steppi would you like @jjerphan to review the math, too, or are you comfortable merging?

The math looks good to me. I’m comfortable with merging.