On Fr, 2014-08-15 at 16:42 +0200, Eelco Hoogendoorn wrote:
> Agreed; this addition occurred to me as well. Note that the
> implemenatation should be straightforward: just allocate an enlarged
> array, use some striding logic to construct the relevant view, and let
> einsums internals act on the view. hopefully, you wont even have to
> touch the guts of einsum at the C level, because id say that isn't for
> the faint of heart...
>
I am not sure if einsum isn't pure C :). But even if, it should be doing
something identical already for duplicate indices on the inputs...
- Sebastian
>
> On Fri, Aug 15, 2014 at 3:53 PM, Sebastian Berg
> <sebast...@sipsolutions.net> wrote:
> On Do, 2014-08-14 at 12:42 -0700, Stephan Hoyer wrote:
> > I think this would be very nice addition.
> >
> >
> > On Thu, Aug 14, 2014 at 12:21 PM, Benjamin Root
> <ben.r...@ou.edu>
> > wrote:
> > You had me at Kronecker delta... :-) +1
> >
>
>
> Sounds good to me. I don't see a reason for not relaxing the
> restriction, unless there is some technical issue, but I doubt
> that.
>
> - Sebastian
>
> >
> >
> > On Thu, Aug 14, 2014 at 3:07 PM, Pierre-Andre Noel
> > <noel.pierre.an...@gmail.com> wrote:
> > (I created issue 4965 earlier today on this
> topic, and
> > I have been
> > advised to email to this mailing list to
> discuss
> > whether it is a good
> > idea or not. I include my original post
> as-is,
> > followed by additional
> > comments.)
> >
> > I think that the following new feature would
> make
> > `numpy.einsum` even
> > more powerful/useful/awesome than it already
> is.
> > Moreover, the change
> > should not interfere with existing code, it
> would
> > preserve the
> > "minimalistic" spirit of `numpy.einsum`, and
> the new
> > functionality would
> > integrate in a seamless/intuitive manner for
> the
> > users.
> >
> > In short, the new feature would allow for
> repeated
> > subscripts to appear
> > in the "output" part of the `subscripts`
> parameter
> > (i.e., on the
> > right-hand side of `->`). The corresponding
> dimensions
> > in the resulting
> > `ndarray` would only be filled along their
> diagonal,
> > leaving the off
> > diagonal entries to the default value for
> this `dtype`
> > (typically zero).
> > Note that the current behavior is to raise
> an
> > exception when repeated
> > output subscripts are being used.
> >
> > This is simplest to describe with an example
> involving
> > the dual behavior
> > of `numpy.diag`.
> >
> > ```python
> > # Extracting the diagonal of a 2-D array.
> > A = arange(16).reshape(4,4)
> > print(diag(A)) # Output: [ 0 5 10 15 ]
> > print(einsum('ii->i', A)) # Same as previous
> line
> > (current behavior).
> >
> > # Constructing a diagonal 2-D array.
> > v = arange(4)
> > print(diag(v)) # Output: [[0 0 0 0] [0 1 0
> 0] [0 0 2
> > 0] [0 0 0 3]]
> > print(einsum('i->ii', v)) # New behavior
> would be same
> > as previous line.
> > # The current behavior of the previous line
> is to
> > raise an exception.
> > ```
> >
> > By opposition to `numpy.diag`, the approach
> > generalizes to higher
> > dimensions: `einsum('iii->i', A)` extracts
> the
> > diagonal of a 3-D array,
> > and `einsum('i->iii', v)` would build a
> diagonal 3-D
> > array.
> >
> > The proposed behavior really starts to shine
> in more
> > intricate cases.
> >
> > ```python
> > # Dummy values, these should be
> probabilities to make
> > sense below.
> > P_w_ab = arange(24).reshape(3,2,4)
> > P_y_wxab = arange(144).reshape(3,3,2,2,4)
> >
> > # With the proposed behavior, the following
> two lines
> > should be equivalent.
> > P_xyz_ab = einsum('wab,xa,ywxab,zy->xyzab',
> P_w_ab,
> > eye(2), P_y_wxab,
> > eye(3))
> > also_P_xyz_ab = einsum('wab,ywaab->ayyab',
> P_w_ab,
> > P_y_wxab)
> > ```
> >
> > If this is not convincing enough, replace
> `eye(2)` by
> > `eye(P_w_ab.shape[1])` and replace `eye(3)`
> by
> > `eye(P_y_wxab.shape[0])`,
> > then imagine more dimensions and repeated
> indices...
> > The new notation
> > would allow for crisper codes and reduce the
> > opportunities for dumb
> > mistakes.
> >
> > For those who wonder, the above computation
> amounts to
> > $P(X=x,Y=y,Z=z|A=a,B=b) = \sum_w P(W=w|
> A=a,B=b) P(X=x|
> > A=a)
> > P(Y=y|W=w,X=x,A=a,B=b) P(Z=z|Y=y)$ with
> $P(X=x|A=a)=
> > \delta_{xa}$ and
> > $P(Z=z|Y=y)=\delta_{zy}$ (using LaTeX
> notation, and
> > $\delta_{ij}$ is
> > [Kronecker's
> >
> delta](http://en.wikipedia.org/wiki/Kronecker_delta)).
> >
> > (End of original post.)
> >
> > I have been told by @jaimefrio that "The
> best way of
> > getting a new
> > feature into numpy is putting it in
> yourself." Hence,
> > if discussions
> > here do reveal that this is a good idea,
> then I may
> > give a try at coding
> > it myself. However, I currently know nothing
> of the
> > inner workings of
> > numpy/ndarray/einsum, and I have higher
> priorities
> > right now. This means
> > that it could take a long while before I
> contribute
> > any code, if I ever
> > do. Hence, if anyone feels like doing it,
> feel free to
> > do so!
> >
> > Also, I am aware that storing a lot of zeros
> in an
> > `ndarray` may not, a
> > priori, be a desirable avenue. However,
> there are
> > times where you have
> > to do it: think of `numpy.eye` as an
> example. In my
> > case of application,
> > I use such diagonal structures in the
> initialization
> > of an `ndarray`
> > which is later updated through an iterative
> process.
> > After these
> > iterations, most of the zeros will be gone.
> Do other
> > people see a use
> > for such capabilities?
> >
> > Thank you for your time and have a nice day.
> >
> > Sincerely,
> >
> > Pierre-André Noël
> >
> _______________________________________________
> > NumPy-Discussion mailing list
> > NumPy-Discussion@scipy.org
> >
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
> >
> >
> >
> > _______________________________________________
> > NumPy-Discussion mailing list
> > NumPy-Discussion@scipy.org
> >
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
> >
> >
> >
> > _______________________________________________
> > NumPy-Discussion mailing list
> > NumPy-Discussion@scipy.org
> > http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
>
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion@scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
>
>
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion@scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
_______________________________________________
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion
