Broadcasting across columns of a matrix


#1

I am struggling to use the dot operator/ broadcasting correctly though when I’d like to broadcast across columns of a matrix.

Suppose I have a function:

function myfun(col_vec, scalar)
   #does something complex operation 
   ...
   return output
end

I’d like a ``vector" version that operates on a matrix by column:

function myfun_vec(mat, vec)
   out = length(output)
   for ix = 1:length(vec)
      out[ix] = myfun(mat[:, ix], vec[ix])
   end
   out
end

Can I use dot syntax / broadcasting to avoid implementing myfun_vec instead? Note

   myfun. ( mat, vec)  

does not return what I want. I mostly want the “.” version so that I can write things like:

myfun. (mat, 2)  #should return myfun_vec(mat, 2 * ones(dim(mat, 2) )

If redefine the data structure of mat to be a vector of vectors and then everything works fine, but that seems super cludgey and doesn’t match surrounding code.


#2

Here is my try,

julia> myfun(col_vec, scalar) = scalar * col_vec
myfun (generic function with 1 method)

julia> mat = collect(reshape(1:10, 2, 5))
2×5 Array{Int64,2}:
 1  3  5  7   9
 2  4  6  8  10

julia> hcat(myfun.([mat[:, i] for i = 1 : size(mat, 2)], collect(1:size(mat, 2)))...)
2×5 Array{Int64,2}:
 1  6  15  28  45
 2  8  18  32  50

I think that in order to use the syntax, myfun.(arg1, arg2), the collections arg1 and arg2 must have the same length.


#3

One option if you know the number of rows of your matrix is to reinterpret it as a vector of SVectors from StaticArrays.jl. This doesn’t copy any data and works very nicely if you have a small number of rows:

julia> A = rand(2, 5)
2×5 Array{Float64,2}:
 0.969289  0.743624  0.280098  0.401384  0.0525567
 0.765344  0.951552  0.913606  0.208708  0.136869 

julia> using StaticArrays

julia> v = reinterpret(SVector{2, Float64}, A);

julia> println.(v)
[0.969289, 0.765344]
[0.743624, 0.951552]
[0.280098, 0.913606]
[0.401384, 0.208708]
[0.0525567, 0.136869]

#4

You can try:

myfun.( [view(mat,:,i) for i in 1:size(mat,2)], vec)

Still, views are allocating but not too much.


#6

Here are some alternatives using mapslices. Both assume that you have a function myfun(col_vec, scalar) as above. Note ii) adds an additional method to myfun().

Cases:
i) The scalars in vec are the same:

mapslices(x->myfun(x,2), mat, dims=1)

ii) The scalars in vec are the same or different e.g. vec = reshape(rand(3), 1,3)

myfun(col_vec) = myfun(col_vec[2:end], col_vec[1])
mapslices(myfun, vcat(vec,mat), dims=1)