kronecker in Stan

[Andrew Gelman]

Hi all. I was speaking with Alex Volfovsky and now I’m all excited again about his hierarchical anova model with Peter Hoff; see here;
http://andrewgelman.com/2013/04/03/hierarchical-array-priors-for-anova-decompositions/
The model uses matrix-normal and array-normal MVN distributions in which the cov matrix is set up as a kronecker product. I recall we had some earlier convos about using kron product in Stan. Anyone interested in this? A possible intermediate solution would be to write a Stan function that implements kron, and if it works we could talk about hardcoding it? What do people think?
See you
A

[Daniel Lee]

I had talked to Peter Hoff about this when I visited UW. He mentioned he had gradients, which means we can implement it in Stan pretty easily! (Not necessarily quickly, but easy.)

I don't have the time to run with it now, but maybe someone else does?

Daniel

A


--
You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Bob Carpenter]

There's a current feature request for Kronecker quadratic forms:

https://github.com/stan-dev/stan/issues/1196

I believe that was to support Gaussian process models.

Andrew --- what do we need to add to Stan in the way of functions
to do what you want to do? If you let me know what we need in
addition to the quadratic form above, I'll create an issue. Anyone
with a GitHub ID can create an issue by clicking on "new issue", so
feel free to do it yourself.

We should talk at the meeting next week about priorities for the
next three to six months.

- Bob

> --
> You received this message because you are subscribed to the Google Groups "stan development mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-dev+u...@googlegroups.com.

[Aki Vehtari]

In that Issue
Bob. "Need utilities for Kroneckers for Gaussian processes."

and

bgoodri:
"I think the most immediate need is to have a user-facing function called something like quad_form_kronecker ... Z' * kronecker_product(X, Y) * Z ..."

I think quad_form_kronecker is useful but not enough.
1) For example for many Gaussian process models we have discussed with Andrew, it would be useful to define a list of small covariance matrices K1, K2, K3, ..., KM instead of just X and Y
2) For GPs, instead of quad form (which you could have, too), it would be better to have Cholesky form which can used with (times vector|left division|log trace|multi_normal_cholesky)
3) With MML it would also useful to be able to just define Kronecker product by just listing the small matrices and never forming the full natrix

Aki

[Ben Goodrich]

We have a matrix normal (log) PDF in C++ already

https://raw.githubusercontent.com/stan-dev/stan/develop/src/stan/prob/distributions/multivariate/continuous/matrix_normal.hpp

but it appears not to be exposed to the Stan language for some reason.

Ben

[Andrew Gelman]

Bob:
I’ll defer to the others on what exactly is needed here. If it can be done initially via a Stan function, then it’s possible that we’re already there, or close to it. One meta-question is: to what extent can we use Stan functions as intermediate septs between “having an idea” and “having it in Stan.”
A

[Seth Flaxman]

Hi Stan. I'm very interested in this; maybe I can collect together the various requests and come to a future Stan dev meeting to help figure out priorities. (Not this coming week, but later is fine.) 

A related question that I expect there'd be much interest in is sparse matrix operations. Is there an issue page for that on the repository / have there already been discussions?

To answer the meta-question: I think with functions Stan can already handle the calculation of a Kronecker matrix-vector product (K1 \otimes K2)v but not (K1 \otimes K2 \otimes K3) v.

Seth

[Bob Carpenter]

[off users list]

We should be having this discussion on stab-dev.  Seth — would you like

me to give you permissions to post there (assuming you are already

subscribed).

Andrew:   We have no mechanism by which to distribute something written

as a Stan function.  We need to figure out includes first.  

The big bottleneck to getting things into Stan is the testing and doc.

Seth  1:  Yes, there’s issues in sparse structures.  There’s better support

in Eigen now.  It’s a pretty major overhaul as it’ll need a design for how to

declare such a structure, then changes to the parser, abstract syntax

tree, code generator, and all the I/O.  Then it’ll percolate through R and

Python and MATLAB and Julia.

Seth 2:  It’d be great if you help us design what we’ll need.

- Bob

You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.

[Aki Vehtari]

It would be nice to have a Kronecker and sparse matrix data types. Then the rest of the Stan code for some model would be same, for example, if the type of covariance matrix is changed (naturally the code constructing the matrix would be different for a Kronecker product matrix).

First you could start by defining the Kronecker data type, but inside functions just compute the corresponding full matrix and then use the normal function. Then you could start modifying functions which can exploit the Kronecker form.

An alternative approach is to have just couple helper functions in Stan language to form Kronecker products and compute quadratic forms etc.

Aki

[Andrew Gelman]

Bob:
Could we have a place on the Github site for contributed functions (where the “contributors” could include people such as you or me or Aki, i.e., “contributed” just means that it hasn’t been fully tested and doc’d)?

[Bob Carpenter]

Sure. We could put it in example-models or we could
create a new public repo called something like stan-functions.
Let’s decide on naming and repo organization and vetting process
for contribs (if any) at the meeting tomorrow.

The drawback for users is that they have to cut-and-paste
them into our models until we get includes sorted out. Then
we could distribute a library of such functions.

The drawback for developers is that we don’t have a good way
of testing Stan functions.

Both of these issues are fairly high priorities for us.

- Bob

[Seth Flaxman]

Is this planned for discussing at Tuesday's meeting? Aki, can you join?

[Aki Vehtari]

I can join depending on time. I'm in Denmark this week.

Aki

[Ben Goodrich]

Did we actually end up switching the meetings to Tuesday mornings EST?

[Andrew Gelman]

I thought that morning switch was for the future, when Bob’s in Australia?