Adding/subtracting a vector from each column of matrix

suppose I have a matrix A and a (column) vector V that has the same height.

Now, if I want to subtract the vector from each column of the matrix, I have two options:

  1. Use a for loop over the columns and do it on each of them
  2. use A - rep_matrix(V, n_cols)

now approach 1 seems the most natural in most languages, but in STAN it apparently creates a lot of unnecessary intermediate values a the autodiff gets rekt

apporach 2 is vectorized and the autodiff seems not so suffer, however it wastes memory and time on allocating the other matrix that isn’t necessary in the cycle approach

So I am thinking that we should perhaps have the following signatures for the arithmetic operators +, -, .*

(matrix, vector) => matrix
(matrix, row_vector) => matrix

which would do the for loop natively avoiding the autodiff problem.

but maybe there is a secret third think the I don’t know about?

1 Like

You’re right that would be optimal—we’d add the analytic derivatives under the hood and simplify the autodiff overhead. We could do that. The place to make feature requests like this is in the issue tracker for the stan-dev/math repository on GitHub (new functions first get implemented there, then imported into the language).

This isn’t quite what we call vectorization (which is itself confusing because it’s not the same vectorization notion as used by CPU optimized code like avx and sse). What’s going on is that it uses matrix addition, which calculates the autodiff all at once. That results in fewer virtual function calls and hence less pointer chasing.

There’s also overhead from looping to copy into the matrix.

We thought about this kind of operation before, but decided that it’d be confusing to just overload our basic operators. So we’d probably add a new named function if we did this.

The best way to implement this is as an expression template in C++ following pattern (2). That will have the effect of copying without actual copying.