MvNormalWishart variate form

[Nimrais]

There are several solutions how the MvNormalWishart issue can be resolved, but let's document our decision somewhere and after I can implement it in a separate PR.

Originally posted by @Nimrais in #147 (comment)

[bvdmitri]

I'm thinking towards something like

abstract type Jointvariate{T} end 

and then

struct MvNormalWishart <: Jointvariate{Tuple{Multivariate, Matrixvariate}}
# ...
end

or something similar.

[Nimrais]

If I understand correctly, the sample is a tuple of a vector and a matrix because this form is more efficient for evaluating the sufficient statistics, right? I am referring to this function getsufficientstatistics(::Type{MvNormalWishart}).

[bvdmitri]

because this form is more efficient

I think it is very reasonable to have it as a tuple (or tuple-like structure of our own), but the entire interface of joint distributions has not been properly reviewed or discussed. Any decision, which is currently in the code, has no real background.

[ismailsenoz]

I have tried to address this issue in the PR #173 by creating a multi-matrix variate type accepting tuple. It is not as generic as Dmitry's proposal but I had to come up with something to test the gradient of the logpartition.

[Nimrais]

I discovered a useful tool for implementing these features at https://github.com/SciML/RecursiveArrayTools.jl. Therefore, we don't need to sample Tuples, and we also won't lose the interfaces for parts. Essentially, these distributions will sample these ArrayPartitions. If we need to access a vector or matrix part of MvNormalWishart, for instance, we can do so through the ArrayPartitions interface, but normally it will behave like a huge vector.