Оптимизировать суммирование Эйнштейна с помощью повторяющихся матриц, но с другим индексом - PullRequest
0 голосов
/ 09 января 2019

Отказ от ответственности: Я сделал хотя о публикации в математический обмен стека или какой-либо другой. Но я слышал от моих друзей по математике, что они не очень часто используют суммирование по Эйнштейну, но я знаю, что машинное обучение часто использует это. Поэтому я разместил эту проблему здесь (как оптимизация производительности алгоритма).

При проведении исследования о матричных вычислениях (например, сколько поэлементных умножений по крайней мере необходимо), я пытался вычислить следующий градиент:

high level equation

, где ABC означает сокращение трех матриц по первой оси (например, 2x3, 2x4 и 2x5 становится 3x4x5 с суммированной осью 2). По сути, он вычисляет градиент нормы 3-матричного сжатия ABC относительно A. Затем снова вычисляет норму этого градиента относительно A.

Это эквивалентно:

traditional writing

или немного упростить (доказано autograd):

enter image description here

Мы можем записать это в форме суммирования Эйнштейна (используется функцией einsum многих пакетов, таких как numpy, tensorflow и т. Д.)

np.einsum('ia,ib,ic,jb,jc,jd,je,kd,ke->ka', A, B, C, B, C, B, C, B, C)

При написании этого я обнаружил, что матрицы B и C повторяются снова и снова при суммировании. Мне интересно, могу ли я упростить эти "много B и C" до некоторой мощности матрицы ? Это должно сделать вещи логарифмически быстрее. Я попытался упростить вручную, но не понял.

Большое спасибо! Пожалуйста, поправьте меня, если что-то, что я говорю, неверно.

1 Ответ

0 голосов
/ 22 января 2019

Я обнаружил, что могу сначала передать (то есть, внешнее произведение) B (формы ib) и C (формы ic), чтобы получить BC тензор формы ibc:

BC = np.einsum('ib,ic->ibc', B, C)

Затем я могу сжать его с помощью транспонирования, чтобы получить квадратную матрицу:

JUMP = np.einsum('ibc,cbj->ij', BC, BC.T)

Эта матрица JUMP является квадратной и действует как оператор «градиентного скачка», который переходит в порядке градиента. Например, если я хочу получить градиент второго порядка, как упомянуто в вопросе, я могу просто умножить это JUMP на градиент первого порядка, чтобы получить это:

g2 = JUMP @ g1

Чтобы получить третий порядок, умножьте квадрат матрицы на него:

g3 = JUMP @ JUMP @ g1

Чтобы получить четвертый ордер, умножьте его на четвертую степень (не куб):

g4 = np.linalg.matrix_power(JUMP, 4) @ g1
...