# JSoC'19: Logdet via *Chebyhutch*

*Chebyhutch*

## Log-determinant

Logarithm determinant, aka *log-det, is *useful in numerous statistical applications and related fields, including machine learning, for example, it is used in calculating log probability density function (logpdf) and maximum log-likelihood function.

The log-determinant of a positive definite matrix

We have previously seen ways to estimate trace of logarithm by using Stochastic Lanczos Quadrature and Chebyshev-Hutchinson Polynomial Interpolation (or Chebyhutch). In this article, we are going to focus on using Chebyhutch to calculate logdet for sparse matrices.

## Chebyhutch for log-determinant

Here is a matrix that I grabbed from SuiteSparse Matrix Collection website. We will be using this for our testing purpose.

Kuu matrix is a Symmetric Positive Definite matrix of size 7102, nonzero count of 340,200 and condition number 1.575800e+04.

using LinearAlgebra, SparseArrays, BenchmarkTools, TraceEstimation, FileIO, JLD2, Distributions, IterativeSolvers

Kuu.jld2A

That's all we need to start our analysis. Let's first take a look at the available method logdet in the LinearAlgebra library.

M = Matrix(A) logdet(M)

Now, we will try running Chebyhutch on the same matrix for log-determinant. Chebyhutch requires 4 parameters:

- A Symmetric Positive Definite Matrix.
- Iteration number (higher iteration number improves trace precision).
- Polynomial Degree (higher degree improves trace accuracy).
- An analytical function such as log or sqrt.

chebyhutch(A, 4, 6, fn = log)

So, it is evident that Chebyhutch can provide a significant boost in the calculation of sparse log-determinant but what about the accuracy? Here is a simple relative error analysis.

act = logdet(M) obv = chebyhutch(A, 4, 6, fn = log) abs((act - obv)/act)

A relative error in the magnitude of 10

## Application - SparseMvNormal

I have already setup the environment for SparseMvNormal, which is internally using Chebyhutch for calculating the logdet. We are going to see a simple application where we find *logpdf *of a sparse covariance matrix with a random mean vector.

# Start by importing the code from environment include("/SparseMvNormalCode") # Generating the vectors for MvNormal example v = rand(size(A, 1)) x = rand(-100.0:2:100.0, size(A, 1)) println("")

Using the above generated v and x, let us first calculate logpdf for MvNormal.

d = MvNormal(v, M) logpdf(d, x)

Now, let us do the same for SparseMvNormal using the original sparse matrix

include("/SparseMvNormalCode") ds = SparseMvNormal(v, A) logpdf(ds, x)

As we can see, the results are good for any numerical purposes.

Since MvNormal uses Cholesky decomposition to calculate logdet we should to do another execution time analysis.

# Wrapper for complete logpdf using MvNormal function clogpdf(A::AbstractMatrix, v::AbstractVector, x::AbstractVector) d = MvNormal(v, A) logpdf(d, x) end # Wrapper for complete logpdf using SparseMvNormal function clogpdf(A::AbstractSparseMatrix, v::AbstractVector, x::AbstractVector) d = SparseMvNormal(v, A) logpdf(d, x) end

clogpdf($M, v, x)

clogpdf($A, v, x)

With this, I will end the SparseMvNormal example.

Thank you for reading!