Умножение массива Numpy вдоль определенной оси
Я хочу умножить два массива numpy разной формы по определенным осям, не меняя их местами вручную и не добавляя «фиктивные» измерения, если это возможно.
Для двух массивов, скажем, A формы (3,2,2) и B формы (3), я хочу получить продукт формы
np.swapaxes(np.swapaxes(A,0,2)*B),0,2)
или альтернативно
A*B[:,None,None]
Однако для больших массивов оба метода, похоже, увеличивают заметное время вычислений по сравнению с начальной инициализацией A формой (2,2,3). Однако я не могу этого сделать, потому что мне нужно, чтобы A имело форму (3,2,2) для последующей операции умножения.
Могу ли я избежать двойной замены осей или добавления размеров «Нет» для сокращения времени расчета?
Обновлено: Привет и спасибо за все ответы, вот более подробный пример того, что я пытаюсь точно вычислить:
import numpy as np
from time import time
A = np.random.rand(100,100,2,20,100)
B = np.random.rand(100)
C = np.random.rand(2,20,100)
D = np.random.rand(100,100,2,20,100)
start = time()
for k in range(98,-1,-1):
A[k,:,:,:,:] = np.pad(C[:,:,k:],((0,0),(0,0),(0,k)),
constant_values=1)*np.pad(B[k:],((0,k)),
constant_values=1)*A[k+1,:,:,:,:] + D[k,:,:,:,:]
end = time()
print('straight forward product: ', end-start)
start = time()
for k in range(98,-1,-1):
A[k,:,:,:,:] = np.pad(B[k:],((0,k)), constant_values=1)
[:,None,None,None]*A[k,:,:,:,:]*np.pad(C[:,:,k:],((0,0),(0,0),
(0,k)), constant_values=1)*A[k+1,:,:,:,:] + D[k,:,:,:,:]
end = time()
print('None Dummy-Dimensions: ', end-start)
второй метод занимает у меня почти вдвое больше времени, но я думаю, что теперь это связано с порядком операторов умножения, а не с тем, какая ось на какую ось умножается.
edit2: Я знаю, что результаты этих двух вычислений не совпадают. В первом случае просто B и C нужно было умножить на A по 5-й оси, поэтому порядок не имел значения, и я случайно выбрал этот. В следующем случае B нужно было умножить на A вдоль 2-й оси, поэтому я изменил порядок умножения, не думая, что это повлияет на производительность, я подозревал, что потеря производительности происходит из-за нотации B[:,None,None,None] .
Можете ли вы предоставить фрагмент кода с правильно определенными A и B?
Ничего плохого в A*B[:,None,None] как таковом нет, это просто стоит мне времени, которое я хотел бы сэкономить на расчетах.
Думаю, можно просто взять np.ones((3,2,2)) для A и np.ones((3)) для B.
Для каких размеров вы видите увеличение времени? B[:,None,None] делает view, что, как показывает обычный опыт, занимает незначительное время. Но размещение очень большого измерения первым, а не последним может привести к разнице во времени из-за проблем с кэшированием/подкачкой памяти.
На практике моя форма A примерно равна (2,20,100,100,100), а B из (100), операцию умножения также пришлось выполнять много раз, например, от 100 до 1000 раз, когда B в моем случае нужно умножить на третью ось.
Ответы 2
В моих быстрых тестах умножение по третьей оси на самом деле происходит быстрее, чем по последней.
In [135]: A=np.ones((2,2,100,100,100)); B=np.ones((100))
In [136]: timeit A*B
17.9 ms ± 333 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
In [137]: timeit A*B[:,None,None]
16 ms ± 1.75 ms per loop (mean ± std. dev. of 7 runs, 100 loops each)
эйнсум
Другой ответ предложил einsum. Я не вижу никаких улучшений
In [138]: timeit np.einsum('ijklm,k->ijklm',A,B)
22.7 ms ± 576 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)
In [139]: timeit np.einsum('ijklm,m->ijklm',A,B)
22.8 ms ± 2.31 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
Спасибо за тестирование. Я как-то недоумеваю, почему нет разницы во времени расчета. Я проверю это еще раз, когда завтра вернусь на работу, чтобы посмотреть, не пропустил ли я что-то еще.
Что касается einsum, я выше отметил, что подозреваю, что это другая операция — сокращение тензора. Имеет ли это смысл, если в вашем тесте потребуется больше времени, например, из-за неявной операции суммирования.
Хотя einsum часто используется для обозначения сокращения, при таком индексировании он ничего sum-of-products не делает. Посмотрите print(np.einsum_path('ijklm,m->ijklm',A,B)[1]). То же самое ijklm 5d как слева, так и справа.
@hpaulj, возможно, вы сможете немного повысить скорость, если при выполнении A*B измените массив на формат Fortran/колонок.
@Меркурий. Это именно то, что проверяет первый тест.
вы можете попробовать использовать функцию библиотеки einsum numpy, которая работает намного быстрее,
import numpy as np
A, B = np.random.rand(3, 2, 2), np.random.rand(3)
%%timeit
np.einsum('ijk,i->ijk', A, B)
>>> 1.19 µs ± 22.3 ns per loop (mean ± std. dev. of 7 runs, 1,000,000 loops each)
На первый взгляд, глядя на функцию einsum, я не уверен, что это то же самое, что и умножение, о котором я упоминал. Разве формы впоследствии не станут другими из-за суммирования? Я ожидаю, что результат вашего примера будет иметь форму (2,2) вместо (3,2,2), не так ли?
@Себастьян, ты путаешь это с "ijk,i->jk". Но это «ijk,i->ijk». Оно останется (3, 2, 2).
Что не так с
A*B[:,None,None]? Я не уверен, что вы можете сделать что-то проще этого. Вы выполняете поэлементное умножение, а не обычное матричное умножение.