Trying to understand broadcast better

Having a list of 2x1 vectors, I wanted to do a for loop for the array and hcat each element with each vector to create a 2x2 matrix. I was able to get what I wanted with a list comprehension but I was wondering if there’s a way to do it with broadcasting hcat.() instread? I guess I’m just used to R and recycling to common vector sizes.

arrays = [[0,0], [1, 1], [1, 2], [2, 1]]
new_array = []
for i in arrays
    transformed = [hcat(i, a) for a in arrays]
    push!(new_array, transformed)
end

How would this work in R? (for those not familiar)
Something like this might work:

julia> arrays = [[0,0], [1, 1], [1, 2], [2, 1]]

julia> [hcat.(Ref(array),arrays) for array in arrays]
4-element Vector{Vector{Matrix{Int64}}}:
 [[0 0; 0 0], [0 1; 0 1], [0 1; 0 2], [0 2; 0 1]]
 [[1 0; 1 0], [1 1; 1 1], [1 1; 1 2], [1 2; 1 1]]
 [[1 0; 2 0], [1 1; 2 1], [1 1; 2 2], [1 2; 2 1]]
 [[2 0; 1 0], [2 1; 1 1], [2 1; 1 2], [2 2; 1 1]]

Y am using Ref() to indicate that broadcasting should not apply to the first argument of the function.

A somehow equivalent way would be using map:

julia> map(array -> hcat.(Ref(array),arrays), arrays)
4-element Vector{Vector{Matrix{Int64}}}:
 [[0 0; 0 0], [0 1; 0 1], [0 1; 0 2], [0 2; 0 1]]
 [[1 0; 1 0], [1 1; 1 1], [1 1; 1 2], [1 2; 1 1]]
 [[1 0; 2 0], [1 1; 2 1], [1 1; 2 2], [1 2; 2 1]]
 [[2 0; 1 0], [2 1; 1 1], [2 1; 1 2], [2 2; 1 1]]
2 Likes

You can use a single comprehension to create the matrix of vectors:

[hcat(a, b) for a in arrays, b in arrays]

But yes you can do it with a simple broadcast: if you broadcast a vector of n elements (or a n \times 1 matrix) with a 1\times m matrix you get an n\times m result:

julia> (1:3) .+ [10 20]
3×2 Matrix{Int64}:
 11  21
 12  22
 13  23

Applying this idea to your problem we can write

hcat.(arrays, permutedims(arrays))

(Note the use of permutedims instead of ' (the latter calculates the adjoint rather than the transpose, and it does it recursively which is not what we want here).

1 Like

Just an aside: Should avoid doing this:

Creating an ‘untyped’ container is a performance killer, and means that your output type is Vector{Any}, which is most likely not what you want.

Every time you see something = [] in your code, it’s a red performance flag.

3 Likes

I just want to add the correct signature:

julia> new_array = Vector{Matrix{eltype(first(arrays))}}[]
Vector{Matrix{Int64}}[]

julia> for i in arrays
           transformed = [hcat(i, a) for a in arrays]
           push!(new_array, transformed)
       end

julia> new_array
4-element Vector{Vector{Matrix{Int64}}}:
 [[0 0; 0 0], [0 1; 0 1], [0 1; 0 2], [0 2; 0 1]]
 [[1 0; 1 0], [1 1; 1 1], [1 1; 1 2], [1 2; 1 1]]
 [[1 0; 2 0], [1 1; 2 1], [1 1; 2 2], [1 2; 2 1]]
 [[2 0; 1 0], [2 1; 1 1], [2 1; 1 2], [2 2; 1 1]]
1 Like

Thanks! I’m trying to learn Julia’s type system more since I’m less experienced in types (coming from Python and R where we tend to skip a lot of that kind of stuff). I’ll be going through the official type documentation.

You can also do something like this, to get the correct type without having to write the signature (and the operation become more generic as well):

julia> for (ind,i) in pairs(arrays)
           transformed = [hcat(i, a) for a in arrays]
           if ind == 1
               new_array = [ transformed ]
           else
               push!(new_array, transformed)
           end
       end

julia> new_array
4-element Vector{Vector{Matrix{Int64}}}:
 [[0 0; 0 0], [0 1; 0 1], [0 1; 0 2], [0 2; 0 1]]
 [[1 0; 1 0], [1 1; 1 1], [1 1; 1 2], [1 2; 1 1]]
 [[1 0; 2 0], [1 1; 2 1], [1 1; 2 2], [1 2; 2 1]]
 [[2 0; 1 0], [2 1; 1 1], [2 1; 1 2], [2 2; 1 1]]

This is because when you initialize an array with one element, it becomes typed for that type of element:

julia> x = []
Any[]

julia> x = [1]
1-element Vector{Int64}:
 1

Or just for fun:

julia> typeof([[[arrays[1];;]]])[]   # Empty array of type given by typeof
Vector{Vector{Matrix{Int64}}}[]
1 Like

Turns out the way I did it in R feels kind of inelegant but :shrug:. Definitely doesn’t look as nice as even just doing a list comprehension and an lapply in an lapply feels weird.

arrays <- list(c(0,0), c(1, 1), c(1, 2), c(2, 1))
result <- lapply(arrays, function(v1) {
  lapply(arrays, function(v2) {
    cbind(v1, v2)
  })
})

Although this code works in the REPL, it won’t work inside of a function:

julia> function make_new_array(arrays)
           for (ind,i) in pairs(arrays)
               transformed = [hcat(i, a) for a in arrays]
               if ind == 1
                   new_array = [ transformed ]
               else
                   push!(new_array, transformed)
               end
           end
           return new_array
       end
make_new_array (generic function with 1 method)

julia> arrays = [[0,0], [1,1], [1,2], [2,1]]
4-element Vector{Vector{Int64}}:
 [0, 0]
 [1, 1]
 [1, 2]
 [2, 1]

julia> make_new_array(arrays)
ERROR: UndefVarError: new_array not defined
Stacktrace:
 [1] make_new_array(arrays::Vector{Vector{Int64}})
   @ Main .\REPL[17]:7
 [2] top-level scope
   @ REPL[19]:1

new_array is not defined on the second and subsequent iterations of the loop. The reason is that loop-local variable new_array is rebound on every loop iteration, as explained in this post.

To fix it, one can add a local declaration:

julia> function make_new_array2(arrays)
           local new_array
           for (ind,i) in pairs(arrays)
               transformed = [hcat(i, a) for a in arrays]
               if ind == 1
                   new_array = [ transformed ]
               else
                   push!(new_array, transformed)
               end
           end
           return new_array
       end
make_new_array2 (generic function with 1 method)

julia> make_new_array2(arrays)
4-element Vector{Vector{Matrix{Int64}}}:
 [[0 0; 0 0], [0 1; 0 1], [0 1; 0 2], [0 2; 0 1]]
 [[1 0; 1 0], [1 1; 1 1], [1 1; 1 2], [1 2; 1 1]]
 [[1 0; 2 0], [1 1; 2 1], [1 1; 2 2], [1 2; 2 1]]
 [[2 0; 1 0], [2 1; 1 1], [2 1; 1 2], [2 2; 1 1]]
4 Likes