diag-plus function

[Andrew Gelman]

Hi all. Sometimes we want to add a diagonal matrix to an existing matrix; for example,
B = A + sigma^2*diag(K)
or maybe
B = A + diag(sigma^2), where sigma here would be a vector of length K.

Aki says that this sort of expression is inefficient because Stan has to waste time including all the +0's in the expression tree; hence instead we add the diagonals in a little loop. This suggests that Stan should have a diag+ function that adds a diagonal vector to a square matrix.

A

[Bob Carpenter]

I created an issue for the primitive diag_add operation, which
would be useful (given that we don't have template expressions yet):

https://github.com/stan-dev/math/issues/411

Looks like it might also be nice to have a scalar version of this, too.

> On Oct 14, 2016, at 11:43 AM, Andrew Gelman <gel...@stat.columbia.edu> wrote:
>
> Hi all. Sometimes we want to add a diagonal matrix to an existing matrix; for example,
> B = A + sigma^2*diag(K)

What is diag(K) if K is an integer? Did you mean something like
unit(K) to create a K x K unit vector? We could use that function, too,
so I added it to the issue.

> or maybe
> B = A + diag(sigma^2), where sigma here would be a vector of length K.

This would be

A + diag_matrix(sigma .* sigma)

in current Stan syntax if sigma is a vector. For now, the thing to do would be
to write a Stan function:

matrix add_diag(vector v, matrix u) {
matrix[rows(u), cols(u)] w;
if (rows(u) != cols(u))
reject("add_diag: matrix must be square");
if (rows(v) != rows(u))
reject("add_diag: vector and matrix must have same num rows");
w = u;
for (k in 1:rows(u))
w[k, k] = v[k];
return w;
}

I threw in the error checking just for Andrew :-) You don't need
it if you are sure u is square and matches the size of v.

- Bob

>
> Aki says that this sort of expression is inefficient because Stan has to waste time including all the +0's in the expression tree; hence instead we add the diagonals in a little loop. This suggests that Stan should have a diag+ function that adds a diagonal vector to a square matrix.
>
> A
>

> --
> 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.
> For more options, visit https://groups.google.com/d/optout.
>
>