Оптимизация операций с матрицей numpy (в настоящее время используется цикл for)

Я написал некоторый код для вычисления n матриц на основе n элементов в списке, а затем в конце умножил все матрицы вместе.

Код относительно медленный, и я хотел бы узнать больше об оптимизации Python. Я использовал инструменты профилирования и определил, что замедление в моей программе связано с циклом умножения матриц.

Интересно, есть ли у кого-нибудь какие-либо рекомендации о том, как я могу ускорить это, возможно, воспользовавшись встроенными функциями на основе C в Python / NumPy?

def my_matrix(x):

    # Initialise overall matrix as an identity matrix
    # | M_11 M_12 |
    # | M_21 M_22 |
    M = np.matrix([[1, 0],[0, 1]])

    for z in z_all:
        param1 = func1(z)
        param2 = func2(x, z)
        param3 = func3(x, z)

        M_11 = param1 + param2
        M_12 = param1 - param2
        M_21 = param1 * param2
        M_22 = param1 / param2

        # Multiply matrix with overall master matrix
        M = M * np.matrix([[M_11, M_12],[M_21, M_22]])
    return M

Из небольшого справочного чтения кажется, что вызовы функций требуют больших вычислительных ресурсов, поэтому вычисление массивов для моих параметров, а затем доступ к массивам может быть более эффективным, чем оценка функции в цикле каждый раз... например.

 param1s = funcs(z_all)
 param2s = funcs(x, z_all)
 etc

а затем в цикле for:

for i, z in enumerate(z_all):
    param1 = params1[i]
    param2 = params2[i]
 etc.

Это быстрее, но только примерно на 10%, поскольку экономия от меньшего количества вызовов функций компенсируется временем, затрачиваемым на доступ к массиву с использованием param1 = params1[i] в ​​цикле.

У кого-нибудь есть рекомендации, пожалуйста?


person IanRoberts    schedule 02.08.2014    source источник
comment
Являются ли какие-либо из ваших параметров/переменных постоянными на протяжении всего цикла? Если это так, вы можете определить его перед циклом и, таким образом, избавить себя от запуска функции снова и снова.   -  person tooty44    schedule 03.08.2014


Ответы (1)


Вы можете векторизовать вычисление M_11, ... M_22, выполнив M_11s = params1 + params2 и т. д.

Таким образом, вам нужно будет выполнять только умножение матриц в цикле:

import numpy as np

...

# compute your 'params' over vectors of z-values
param1s = func1(z_all)
param2s = func2(x, z_all)
param3s = func3(x, z_all)  # you don't seem to be using this for anything...

# compute 'M_11, ... M_22'
M_11 = param1s + param2s
M_12 = param1s - param2s
M_21 = param1s * param2s
M_22 = param1s / param2s

# we construct a (2, 2, nz) array from these
M_all = np.array([[M_11, M_12], [M_21, M_22]])

# roll the 'nz' axis to the front so that its shape is (nz, 2, 2)
M_all = np.rollaxis(M_all, -1, 0)

# initialize output with the identity
M_out = np.eye(2)

# loop over each (2, 2) subarray in 'M_all', update the output with the
# corresponding dot product
for mm in M_all:
    M_out = M_out.dot(mm)
person ali_m    schedule 02.08.2014
comment
Спасибо, это дает ускорение в 3 раза! Не могли бы вы уточнить, почему это быстрее, пожалуйста? Лучше ли работать с массивами, чем с матрицами в целом, поскольку точечный продукт numpy здесь кажется эквивалентным матричному умножению. - person IanRoberts; 03.08.2014
comment
Основная причина ускорения заключается в том, что вычисление M_11, ... ,M_22 теперь векторизовано. Работа с целыми векторами всегда предпочтительнее индексации отдельных элементов. Основной аргумент в пользу использования здесь массивов, а не матриц, заключается в том, что массивы могут иметь произвольный ранг, тогда как максимальный ранг матриц равен 2. Если бы я использовал матрицы, M_all должен был быть списком матриц (2, 2), а не одним массивом. Есть много других причин предпочесть массивы матрицам, см. здесь. - person ali_m; 03.08.2014
comment
Спасибо, это действительно полезный ответ. На самом деле я хочу вычислить M_out для 2D-пространства параметров. На данный момент я перебираю один параметр (назовем его временем) и для каждого времени в массиве времен я использую приведенный выше код для вычисления M_out. Интересно, можно ли это улучшить, заменив цикл for некоторой векторной операцией. Я разработал и обобщил свой вопрос, так как считаю, что слишком много обсуждений здесь может разбавить этот вопрос и ответ: stackoverflow.com/questions/25100046/ - person IanRoberts; 03.08.2014