Fix Gaussian States

[charlesknipp]

For now, linear models are defined using ::PDMat which not only creates extraneous allocations, but also requires a Cholesky decomp upon object construction.

Since sampling from some linear models is not possible by the typical methods defined by Distributions.jl and Gaussian.jl, we have to use this modified Cholesky process for matrices with zeros on the diagonal. The only obvious choice was to make sampling easy by defining a PDMat with the StateSpaceModels.jl's cholesky_decomposition(x::AbstractMatrix).

An alternative approach would be to drop the PDMats dependency and instead construct a new class of distributions ala DistributionsAD.jl, where MvNormals are defined by the upper Cholesky of the covariance, and sampling random normals is done with CuArrays.

If I wanted to go the extra mile, it may be worth using the GPU to do some computations. Although, this isn't a priority

[charlesknipp]

The newest commit demonstrates yet another problem with linear models. Minor changes to harvey-trimbur.jl, which define the "balanced cyclical model", break the default constructor for linear models. This is not unexpected since the constructor was hastily programmed to begin with.

Cholesky fails when trying to construct the PDMat form of the trend covariance. This can be fixed directly by defining the known Cholesky decompositions of the block diagonals.

import StateSpaceInference: cholesky_decomposition

# cycle covariance
Σψ = kron(diagm([zeros(n-1)..., 1]), σκ²*I(2))
Σψ = PDMat(cholesky_decomposition(Σψ))

# trend covariance
Σx = UpperTriangular(sqrt(σε²)*ones(m, m))
Σx = PDMat(Cholesky(Σx))

However, this is a convenience of the model. Ideally, we should take advantage of dispatch to have a one size fits all type solution. Theoretically, this decomposition exists even for the provided matrix types; but as we've seen, this is not always what Julia determines.

• Idea 1: we can instead define linear Gaussian models by using the Cholesky factorization of the covariance.
• Idea 2: we can continue to account for more scenarios within the cholesky_decomposition function.

I much prefer the former, because cholesky_decomposition exists purely to construct PDMats from matrices with zeros in the diagonals. The less I rely on PDMats the better.