numpy - обозначение einsum: скалярное произведение стека матриц со стеком векторов

Question

Я хочу умножить n-dim стек из m * m матриц на n-dim стек векторов (длина m), чтобы результирующий массив m * n содержал результат точечного произведения матрицы и векторав n-й записи:

vec1=np.array([0,0.5,1,0.5]); vec2=np.array([2,0.5,1,0.5])
vec=np.transpose(n.stack((vec1,vec2)))
mat = np.moveaxis(n.array([[[0,1,2,3],[0,1,2,3],[0,1,2,3],[0,1,2,3]],[[-1,2.,0,1.],[0,0,-1,2.],[0,1,-1,2.],[1,0.1,1,1]]]),0,2)
outvec=np.zeros((4,2))
for i in range(2):
    outvec[:,i]=np.dot(mat[:,:,i],vec[:,i])

Вдохновленный этим постом Мудрое произведение элементов на матрицы и векторы , я попробовал все различные возмущения комбинаций индексов в einsum и обнаружил, что

np.einsum('ijk,jk->ik',mat,vec)

дает правильный результат.

К сожалению, я действительно не понимаю этого - я предположил, что повторение записи k в части 'ijk, jk' означает, что я умножаю ANDсумма по к.Я пытался прочитать документацию https://docs.scipy.org/doc/numpy-1.15.1/reference/generated/numpy.einsum.html,, но все еще не понимаю.

(Мои предыдущие попытки включали,

 np.einsum('ijk,il->ik', mat, vec)

Я даже не уверен, чтоэто значит. Что будет с индексом l, когда я его урону?)

Заранее спасибо!

Answer 1

Читайте Обозначение Эйнштейна * .

В основном, применяются следующие правила:

Без ->

• Любая буква, повторяемая во входных данных, представляет собой ось, которая должна быть умножена и суммирована по
• Любая буква, не повторяемая на входах, включается в вывод

С ->

• Любая буква, повторяемая на входах, представляет собой ось, которая должна быть умноженаover
• Любая буква, отсутствующая в выходных данных, представляет собой ось, которая будет суммироваться по

Так, например, с матрицами A и B с одинаковой формой:

np.einsum('ij, ij',       A, B)  # is A ddot B,                returns 0d scalar
np.einsum('ij, jk',       A, B)  # is A dot  B,                returns 2d tensor
np.einsum('ij, kl',       A, B)  # is outer(A, B),             returns 4d tensor
np.einsum('ji, jk, kl',   A, B)  # is A.T @ B @ A,             returns 2d tensor
np.einsum('ij, ij -> ij', A, B)  # is A * B,                   returns 2d tensor
np.einsum('ij, ij -> i' , A, A)  # is norm(A, axis = 1),       returns 1d tensor
np.einsum('ii'             , A)  # is tr(A),                   returns 0d scalar

Answer 2

einsum легко (если вы некоторое время играли с перестановкой индексов, то есть ...).

Давайте поработаем с чем-то простым, тройным стеком 2 × 2 матриц и тройной стек 2 ×, массивов

import numpy as np

a = np.arange(3*2*2).reshape((3,2,2))
b = np.arange(3*2).reshape((3,2))

Нам нужно знать, что мы собираемся вычислить, используя einsum

In [101]: for i in range(3): 
     ...:     print(a[i]@b[i])                                                                            
[1 3]
[23 33]
[77 95]

Что мы сделали?у нас есть индекс i, который фиксируется, когда мы выполняем скалярное произведение между одной из суммированных матриц и одним из суммированных векторов (оба индексированы i), а отдельная строка вывода подразумевает суммирование по последнему индексуматрица с накоплением и одиночный индекс с накопленным вектором.

Это легко кодируется в директиве einsum

• , мы хотим, чтобы тот же индекс i определял матрицу,вектор, а также выходные данные,
• , которые мы хотим уменьшить по последнему индексу матрицы и оставшемуся векторному индексу, скажем, k
• , мы хотим, чтобы на выходе было столько столбцов, сколько строкв каждой суммированной матрице, скажем, j

Следовательно

In [102]: np.einsum('ijk,ik->ij', a, b)                                                                   
Out[102]: 
array([[ 1,  3],
       [23, 33],
       [77, 95]])

Я надеюсь, что мое обсуждение того, как я получил директиву правильно, ясно, правильно и полезно.

Answer 3

In [321]: vec1=np.array([0,0.5,1,0.5]); vec2=np.array([2,0.5,1,0.5])
     ...: vec=np.transpose(np.stack((vec1,vec2)))
In [322]: vec1.shape
Out[322]: (4,)
In [323]: vec.shape
Out[323]: (4, 2)

Приятной особенностью функции stack является то, что мы можем указать ось, пропуская транспонирование:

In [324]: np.stack((vec1,vec2), axis=1).shape
Out[324]: (4, 2)

Почему сочетание np. и n.?NameError: name 'n' is not defined.Такие вещи почти отсылают меня прочь.

In [326]: mat = np.moveaxis(np.array([[[0,1,2,3],[0,1,2,3],[0,1,2,3],[0,1,2,3]],[[-1,2.,0
     ...: ,1.],[0,0,-1,2.],[0,1,-1,2.],[1,0.1,1,1]]]),0,2)
In [327]: mat.shape
Out[327]: (4, 4, 2)

In [328]: outvec=np.zeros((4,2))
     ...: for i in range(2):
     ...:     outvec[:,i]=np.dot(mat[:,:,i],vec[:,i])
     ...:     
In [329]: outvec
Out[329]: 
array([[ 4.  , -0.5 ],
       [ 4.  ,  0.  ],
       [ 4.  ,  0.5 ],
       [ 4.  ,  3.55]])

In [330]: # (4,4,2) (4,2)   'kji,ji->ki'

Из вашей петли расположение оси i (размер 2) понятно - последнее во всех 3 массивах.Это оставляет одну ось для vec, давайте назовем это j.Он соединяется с последним (рядом с i из mat).k переносится с mat на outvec.

In [331]: np.einsum('kji,ji->ki', mat, vec)
Out[331]: 
array([[ 4.  , -0.5 ],
       [ 4.  ,  0.  ],
       [ 4.  ,  0.5 ],
       [ 4.  ,  3.55]])

Часто строка einsum записывается сама.Например, если mat было описано как (m, n, k) и vec как (n, k), с результатом (m, k)

В этом случае только jизмерение суммируется - оно появляется слева, но справа.Последнее измерение, i в моей записи, не суммируется, потому что, если оно появляется с обеих сторон, так же, как и в вашей итерации.Я думаю об этом как о том, что нужно ехать.Он не является частью продукта dot.

По сути, вы накладываетесь на последнее измерение, размер 2.Обычно мы ставим на первое, но вы транспонируете оба, чтобы поместить это последнее.

---

Ваша «неудачная» попытка выполняется и может быть воспроизведена как:

In [332]: np.einsum('ijk,il->ik', mat, vec)
Out[332]: 
array([[12. ,  4. ],
       [ 6. ,  1. ],
       [12. ,  4. ],
       [ 6. ,  3.1]])
In [333]: mat.sum(axis=1)*vec.sum(axis=1)[:,None]
Out[333]: 
array([[12. ,  4. ],
       [ 6. ,  1. ],
       [12. ,  4. ],
       [ 6. ,  3.1]])

j и l размеры не отображаются справа, поэтому они суммируются.Их можно суммировать до умножения, потому что они появляются только в одном слагаемом.Я добавил None для включения трансляции (умножив ik на i).

np.einsum('ik,i->ik', mat.sum(axis=1), vec.sum(axis=1))

---

Если вы сложили на первом, и добавил измерение для vec(2,4,1), это будет matmul с (2,4,4) мат.mat @ vec[...,None].

In [337]: m1 = mat.transpose(2,0,1)
In [338]: m1@v1[...,None]
Out[338]: 
array([[[ 4.  ],
        [ 4.  ],
        [ 4.  ],
        [ 4.  ]],

       [[-0.5 ],
        [ 0.  ],
        [ 0.5 ],
        [ 3.55]]])
In [339]: _.shape
Out[339]: (2, 4, 1)