Unintuitive Julia result: selecting a row from a matrix

I bumped into this example trying to convert some julia to python and was surprised. In the example here, python and matlab give the same answer, while julia differs…

If you start with…
X=rand(3,100) (and do X=mat(X) if in python)

and then do the following…

matlab/octave: X(1,:)*X(1,:)'
julia: X[1,:]*X[1,:]'
python: X[1,:]*X[1,:].T

what size should the result be and why? in python and matlab (and julia), X[1,:] will select a row. python and matlab however return a row vector (or a 1xN matrix) while julia returns a column vector (Nx1). there may be a reason, but I am not sure what that is. Is this expected?

2 Likes
julia> x[1,:] # this returns a Vector, which is a column

julia> x[1,:]' # this traspose the Vector into a row-vector
1 Like

Yes, it’s expected. X[1,:] gives a flat vector. In case you want an 1xn array, do X[1:1,:].

4 Likes

Yup, there absolutely is a reason and it’s a doozy. It turns out that the embedding of one-dimensional containers into the algebra of matrices is not obvious. There are a number of choices you must make in how they behave, impacting everything from indexing to transposes. It’s one of the few design decisions that has had so much thought put into it that many Julia folk who have been around for a while have it memorized (#4774). It’s spawned a paper (still a draft) and presentations. There’s even an :4774: emoji :4774: commemorating its author and major thought-leader here and on Slack.

The end-result, though, is a remarkably simple rule: The resulting dimensionality of an indexing expression is the sum of the dimensionalities of the indices. Or even more specifically, the resulting axes of an indexing expression is the concatenation of all the axes of all the indices used.

So that’s why X[1, :] is a 1-dimensional vector while X[[1], :] is a 1xN row matrix.

30 Likes

I hope this is in a very visible place in docs and I just happened to miss it :expressionless:

2 Likes

It appears here: Multi-dimensional Arrays · The Julia Language, which is where we talk about all things indexing.

3 Likes

no but your summary here can beat 3 paragraphs and lengthy example in the actual manual I think.

5 Likes

One thing to keep in mind, that may make this less weird, is that Julia does not return an Nx1 array, but simply an N vector (no ‘x1’). The distinction makes a difference.

6 Likes

This is actually not quite right. First of all, numpy does not return a row vector, but a flat vector like Julia:

In [21]: X[1,:].shape
Out[21]: (5,)

It prints like a row, but it’s not.

Also, you should note that Matlab and numpy disagree with each other (as well as with Julia) in other ways:

Numpy:

  • ‘row slice’ produces a flat vector
  • transposing a flat vector does nothing, the vector is unchanged
  • multiplication with * is element-by-element and broadcasting, not a linear algebra operation
  • => X[1,:]*X[1,:].T returns an N vector

Matlab

  • ‘row slice’ produces a 1xN matrix, not a vector
  • transposing it creates an Nx1 matrix
  • multiplication with * is a linear algebra operation
  • => X(1,:)*X(1,:)' returns a 1x1 matrix (‘scalar’-like)

Julia:

  • ‘row slice’ produces a flat vector
  • transposing it produces a special transpose (or in the case of ' an adjoint) vector, which is ‘row-like’
  • multiplication with * is a linear algebra operation
  • => X[1,:]*X[1,:]' returns an NxN matrix

So it’s not like Julia is the odd one out here. They all disagree.

7 Likes

I noticed this just now, and it’s not something I’ve seen before. Does it alter the behavior of numpy arrays?

Yes, we need to pester @jiahao to make his draft paper public; if I remember right he details the “universes” of reasonable ways you can get vectors to work with matrices. None are perfect.

4 Likes

very interesting. I can see the reason. I caught me by surprise because of the difference with matlab and python with matrix operations. python uses two data types (array and matrix) so the operations for matrix are a bit different. matlab only had matrix operations. at least now I can understand the logic. thanks!

1 Like

yes it does. it make all operations 2d – so that you’re only doing matrix operations. in this way it behaves more like matlab, which is why it was introduced. you’re correct that for an array X[1,:] returns a flat array.

1 Like

Just be careful here, x[1,:]' will lead to a conjugate transpose in case x is complex.

2 Likes