HighDimensionalEconLab / DifferentiableStateSpaceModels.jl

MIT License
46 stars 1 forks source link

Incomplete rank cholesky with pivoting #136

Closed jlperla closed 2 years ago

jlperla commented 2 years ago

It seems like the pivoted cholesky does the trick for dealing with ergodic distributions that aren't full rank. THis can be done by passing in Val(true) and turning off the check.

The problem is that because it is pivoted it has a permutation matrix and you can't directly use the existing code. Not sure if it is possible to get this inside of PDMats (and hence the MvNormal) with the current PDMats and Distributions code?

Alternatively, the kalman could be written to pass in and work with the pivoted cholesky in https://github.com/SciML/DifferenceEquations.jl/blob/main/src/algorithms/linear.jl but it is a major consideration on the right algorithms/etc.

To see the pivotingg in action:

using LinearAlgebra, PDMats, Distributions, Test
A_rank = 2
N = 3
x = rand(N, A_rank)
A = x * x' # rank = 2
Apiv = cholesky(A, Val(true); check = false)
U = Apiv.U[1:(Apiv.rank), invperm(Apiv.p)]  # even though lower rank, using the full rank
@test A ≈ U' * U

A_chol = Cholesky{Float64,Matrix{Float64}}(zeros(3, 3), 'U', 0)
A_chol.factors[1:(Apiv.rank), :] .= U
@test A ≈ A_chol.L * A_chol.U