Re: [Numpy-discussion] What's wrong with matrices?

[Bill B. <wb...@gm...>]

On 7/7/06, Ed Schofield <sch...@ft...> wrote:
>
> Bill Baxter wrote:
> > I would be all for a matrix class that was on equal footing with array
> > and as easy to use as matrices in Matlab.  But my experience using
> > numpy.matrix was far from that, and, given the lack of enthusiasm for
> > matrices around here, that seems unlikely to change.  However, I'm
> > anxious to see what Ed has up his sleeves in the other thread.
>
> Okay ... <Ed rolls up his sleeves> ... let's make this the thread ;)
> I'd like to know why you, Sven, and anyone else on the list have gone
> back to using arrays after trying matrices.  What was inconvenient about
> them?  I'd like a nice juicy list.  The whole purpose of the matrix
> class is to simplify 2d linear algebra.  Where is it failing?


Okay, here are a few that come to mind.
1) Functions that take a matrix but return an array.  Maybe these are all
fixed now.  But they better be fixed not just in numpy but in scipy too.  To
me this implies there needs to be some standard idiom for how to write a
generic array-protocol-using function so that you don't have to think about
it.

2) At the time I was using matrix, scalar * matrix was broken.  Fixed now,
but that kind of thing just shouldn't happen.  There should be a tests for
basic operations like that if there aren't already.

3) mat() doesn't make sense as a shortcut for matrix construction.  It only
saves 3 letters over typing matrix(), and asmatrix is generally more
useful.  So mat() should be a synonym for asmatrix()

4) eye,empty,rand,ones,zeros,arange and anything else that builds an array
from scratch or from a python list should have a matrix equivalent.

5) I've got squeezes like crazy all over my matrix-using code.  Maybe this
was a bug in 0.9.5 or so that's been fixed?  I do seem to recall some
problem with indexing or c_ or something that was causing matrices to grow
extra levels of length 1 axes.  Again, like the scalar*matrix bug, things
like that shouldn't happen.

6) No good way to do elementwise operations?  Sometimes you just want to do
an elementwise mult or divide or exponentiation.  I guess you're supposed to
do  Z = asmatrix(X.A * Y.A).  Yah, right.

7) Finally, once all that is fixed, I find the slavish adherance to ndim=2
to be too restrictive.
    a) Sometimes it's useful to have a stack of matrices.  Pretty often, in
fact, for me.  I guess I could use a python list of matrices or an object
array of matrix or something, but I also think there are times when it's
useful to treat different pairs of axes as the 'matrix' part.   So I'd like
matrices to be able to have ndim>2.
    b) On the other end, I think ndim<2 is useful sometimes too.  Take a
function like mean(), for example.  With no arguments  the return value is a
1x1 matrix (as opposed to a scalar).  That just looks silly.  And it doesn't
work everywhere a scalar would, e.g. with a python built-in like round().
Or take indexing.  It seems odd to me that where() reurns a tuple of
shape==(1,N) objects instead of just (N,) .  Maybe I can get over that
though, as long as it works for indexing (which it seems it does).   But I
think the scalar return case is a real issue.  Here's another: sum().  For
an array you can do sum(sum(a)) and get a scalar if a is 2-d, but for matrix
sum(sum(m)) is the same as sum(m).   And along these lines, m[newaxis] just
silently doesn't do anything.  That doesn't seem right.

That's all I can think of for now.

Beyond that, looming in the back of my mind there's the question of what
happens when we decide we want to have sparse matrices?  Or masked sparse
matrices, for that matter.  These specialty array objects all seem to me
like they should be some sort of mix-ins.  The features they provide are
(conceptually, anyway) independent.  Concrete sibling classes doesn't really
seem like the right inheritance relationship.


> I'd like to help to make matrices more usable.  Tell me what you want,
> and I'll work on some patches.


Thanks for giving it a shot.  I'm currently trying to compile numpy from svn
on Win32 for the first time myself.  If it goes well I can probably help run
tests or maybe even on some coding.  But I still have yet to get past ATLAS.

--Bill