`eachrow` over Array{Any, 2} does not return rows, why?

General Usage
1

I use Julia often, but mostly JuMP. I have incentivized the use of Julia in my laboratory. Today, one of my laboratory colleagues has said to me he used Julia to make a small program, and while he overall liked the language he found the inconsistency between row and column major to have made him lost some sanity and wasted considerable development time. I was confused because I do not make much use of matrices but this somewhat did not align with I had seen in many discussion about Julia. After 15~30 minutes I found where he got confused, and I kinda conceded the point (but I would like to hear you opinion about it). The problem is that two approaches that seem like they should work seamlessly do not work at all. In the examples below, we want to compute a dot product.

julia> m = [1 2 3; 4 5 6]
2×3 Array{Int64,2}:
 1  2  3
 4  5  6
julia> for row in eachrow(m)
    row * [2, 2, 2]
end

As Vectors/Array{Any, 1} are seen as columns, one would expect that to work, but instead you see an error message starting with:
ERROR: MethodError: no method matching *(::SubArray{Int64,1,Array{Int64,2},Tuple{Int64,Base.Slice{Base.OneTo{Int64}}},true}, ::Array{Int64,1})

What is kinda reasonable but then it is followed by a ‘closest candidate’ for the multiplication operator that has two parameters each with literal 20 lines of type description (::LinearAlgebra.Adjoint{#s627,#s626} where #s626<:Union{DenseArray{T,2}, Base.ReinterpretArray{T,2,S,A} where S ...).

If you instead do:

julia> col_factor = fill(2, 3, 1)
3×1 Array{Int64,2}:
 2
 2
 2

julia> for row in eachrow(m)
    row * col_factor
end

You get the much more legible error message saying: ERROR: DimensionMismatch("matrix A has dimensions (3,1), matrix B has dimensions (3,1)"). However, this is also confusing because, in this case, the programmer has no reason to believe that eachrow does not actually returns rows but columns, they may fist guess that they inverted the order of the dimensions in the fill call or that the error is bubbling up in a place different from where the real mismatch happens.

So, was this intended? an oversight? or seems good in theory but there other more general issues that prevent changing this behavior? or I am missing something and one should not expect what my friend expected? Seems to me that it could be a subtle source of bugs.

To clarify even more, why:

julia> for row in eachrow(m); println(row); end
[1, 2, 3]
[4, 5, 6]

and not

julia> for row in eachrow(m); println(row); end
[1 2 3]
[4 5 6]
2

About the dot product, please see my comment here

The rest of this post is about broadcasting and using eachrow, eachcol:

You wrote row * [2, 2, 2] where you meant row .* [2, 2, 2].
To multiply two scalars, use *. To multiply a vector use .*, it is called broadcasting and is a very powerful means of expression that Julia offers.
see dot-syntax and then broadcasting.

So, these work:

julia> m = [1 2 3; 4 5 6]
2×3 Array{Int64,2}:
 1  2  3
 4  5  6

julia> for row in eachrow(m)
          println(row .* [2, 2, 2])
       end
[2, 4, 6]
[8, 10, 12]

julia> for row in eachrow(m)
           println(row .* 2)
       end
[2, 4, 6]
[8, 10, 12]

julia> for col in eachcol(m)
           println(col .* 3)
         end
[3, 12]
[6, 15]
[9, 18]

julia> for col in eachcol(m)
           println(col .* [2,3])
         end
[2, 12]
[4, 15]
[6, 18]

eachrow and eachcol are iterators, to get their content directly, use collect with “broadcasting”:

julia> collect.(eachrow(m)) # notice the `.`
2-element Array{Array{Int64,1},1}:
 [1, 2, 3]
 [4, 5, 6]

julia> collect.(eachrow(m) .* 2)
2-element Array{Array{Int64,1},1}:
 [2, 4, 6]
 [8, 10, 12]
1 Like
3

…

I think you misunderstood my problem. I am sure I did not meant row .* [2, 2, 2], I know what broadcasting is, and is not what I want here. I said I wanted the dot product, also called inner product sometimes. For an example:

julia> [1 1 1] * [1, 1, 1]
1-element Array{Int64,1}:
 3
julia> [1 2 3] * [1, 2, 3]
1-element Array{Int64,1}:
 14
4

If I understand correctly, your question is why eachrow returns a 1x3 and not a 3x1, correct?

1 Like
5

My intuition is actually that the each row of a matrix should be a vector but not 1 times n row, (and each column of a matrix should be a vector too but not a n*1 matrix) which is the same with the behavior here. And this behavior is quite consistent, for example,

julia>  m = [1 2 3; 4 5 6]
2×3 Array{Int64,2}:
 1  2  3
 4  5  6

julia> m[1,:]
3-element Array{Int64,1}:
 1
 2
 3

And each row is indeed a vector

julia> for row in eachrow(m); println(size(row)); end
(3,)
(3,)

And I also expect that julia should give a matrix when I want to get multiple rows together, so I also try the following and get the expected result:

julia> m[1:1,:]
1×3 Array{Int64,2}:
 1  2  3
1 Like
6

to use the dot product:

julia> using LinearAlgebra
julia> m = [1 2 3; 4 5 6];
julia> for row in eachrow(m)
           println( dot(row, [2, 2, 2]) )
       end
12
40

The dot product operator is called dot and it lives in the LinearAlgebra standard library.

1 Like
7

I’ll admit this is actually a bit confusing to me as well. This is actually one reason why I just use dot any time I am multiplying vectors (or row vectors, basically something with a singleton dimension) – no need to care about direction then!

1 Like
8

I think a good way to make this less confusing would be to have a better error message when multiplying two AbstractVectors. It could point users in the right direction by mentioning broadcasted multiply, dot and maybe adjoint/transpose.

3 Likes
9

1d slices in Julia always return a 1d array (interpreted as a column vector), whether you are slicing rows or columns. That is, matrix[:, j] and matrix[i, :] both return 1d arrays. eachrow is consistent with this. In general, any singleton dimension is dropped in a slice operation.

If you want the dot product of two column vectors (1d arrays), do dot(x, y) or x' * y or x'y.

7 Likes
10

@simeonschaub This would be nice. It is just a question of defining Base.:*(a :: AbstractVector, b :: AbstractVector) with a body that just gives a better error message, no? The question would be if it should live in the standard library, or if it should be part of a package that enhances the experience for newcomers.

@tbeason @JeffreySarnoff I was not aware of the dot operator as I do not do much linear algebra, it is nice that LinearAlgebra is in the standard library, if it was not could be kind overkill to install a package just for this.

@Non-Contradiction I kinda agree that in many cases dropping dimensions is probably the most sane behavior. In a construct with four or five dimensions I probably would not expect/want/like to get a 1x1x3x1x1 subarray; but I would understand why I did get it, and would be probably trivial to drops all unitary dimensions (i.e., there would exist a standard library method for this). The point that I agree with my friend is that: if there is a pair of methods called eachrow and eachcol that work over two-dimensional matrices, I would expect without even explicitly thinking about it that eachrow return rows, and not one-dimensional-arrays-of-a-column-major-language which is the same as the eachcol returns. I am happy that the behavior is consistent, but I think the methods will probably lure unaware programmers and cause headaches to them, I think that maybe they should not exist or then be more explicit about this behavior (the documentation just say it returns views). I also think I would not bat an eye for this behavior coming from eachslice but it is the explicit reference to rows/columns that cause confusion.

@stevengj I think this is sane for slicing. However, in methods called eachrow and eachcol maybe the unitary dimension dropping should need to be explicited.

All this said, my only interest in questioning this decision is to understand its reason, and if it will/should stay this way. The idea is that I may actually argue why this is a good design if I find another newcomer with this problem, or that this is changed and no newcomers even have this problem again. The pointer to dot is already something to give as an ‘workaround’.

11

There has been a ton of discussion and design work over a number of years to get here with arrays. Any alteration, even of a minor sort cannot be introduced prior to Julia 2, given the assurances about the stability of Julia 1.

12

I did not suggest: (1) That any breaking change should be done before version 2.0; (2) that would be a change in the design of arrays.

Also, for an example of something I suggested, removing eachrow and eachcol would not be so big of a problem as they were added in Julia 1.1 (see Julia v1.1 Release Notes ¡ The Julia Language) and therefore they are not part of the long term Julia 1.0 version.

13

Thanks, I did not mean to suggest your perspective, merely to inform others.

14

Just to clarify the versioning point— if they’re added in 1.1 and part of the public API (as eachrow and eachcol are) they can’t be removed until 2.0. That’s because any update between 1.1 and 2.0 isn’t a major version change and so shouldn’t be breaking. So it’s not just what’s in 1.0 that matters.

3 Likes
15

You are correct. I expressed myself badly. What I meant was that there
would never be a long term version supporting these methods if they
decided to remove them at Julia 2.0. So, concerning the people only
using the LTS versions, they would never touch them.

1 Like
16

I am not sure what that error message would be. Do you have a suggestion?

I think that the current error message telling you there is no such method is perfectly fine. Generally, it is rarely worthwhile to do anything else, because it would require anticipating user expectations which is usually impossible or difficult. Deprecations and similar are of course an exception.

17

Personally, I think it would be great if it returned something like

ERROR: MethodError: no method matching *(::Array{Int64,1}, ::Array{Int64,1}).
Did you mean `LinearAlgebra.dot(a, b)` or `a'*b` (to obtain a scalar), `a.*b` (for elementwise multiplication), or `a*b'` (for an outer product)?
Closest candidates are:
...

For anyone looking to get her/his feet wet with a Julia PR, this would be an excellent first issue! There’s even a recent example you can copy. The initial issue report is

https://github.com/JuliaLang/julia/issues/34636

and the fix was submitted in

https://github.com/JuliaLang/julia/pull/34642

You could mimic that fix very directly.

8 Likes
18

I understand that each guess about what the user intended to do would help someone, but at the same time it is a small increase in the maintenance burden (all code is).

I think it would be better to figure out a rough general guideline about what to include first. Eg should b * A, eg

ones(2) * ones(2, 2)

error with suggestions for A * b and b' * A?

19

@tim.holy Great! My first PR was something similar. I will check with my lab colleagues if someone wanna do it, if not, I can do it.

@Tamas_Papp In this case, I really think the extra maintenance cost (which is minimal for a method that just printing an error message) is a small price to pay to all hours of the newcomers that find this bug, and of regulars that possibly help those newcomers. A guideline would be great, but it probably should be assembled by the standard library maintainers. For now, they will probably approve or not PRs in a one-by-one case.

2 Likes
20

The cost is not one-time: these messages have to be kept up to date, preserved when the code changes (eg in case of a refactoring), unit tested, and maintained.

If we do this consistently, all kinds of suggestions from the language can be reasonable requested whenever someone is surprised about something.

I am not sure why you think this is a bug, it isn’t — as pointed out above by @stevengj, it is consistent with how slicing works in Julia.

I am not trying to be a curmudgeon here, I just think that it is unreasonable to expect no surprises when learning a new language. People come with prior experience in a whole lot of other languages with conflicting syntax and semantics, so it is very difficult to guess what someone intended in a specific context (eg see #34706 recently, there are many other examples).

These are tiny things that people get over very quickly once they adapt to Julia, and effort is better spent on adding examples to docstrings, and making interfaces or consistent.

2 Likes