Произведение последовательности в NumPy


Мне нужно реализовать следующую функцию с помощью NumPy -

Введите описание изображения здесь

Где F_l(x) - это N число массивов, которые мне нужно вычислить, которые зависят от массива G(x), который мне задан, и A_j - это N коэффициенты, которые также заданы. Я хотел бы реализовать его в NumPy, так как мне пришлось бы вычислять F_l(x) для каждой итерации моей программы. Фиктивный способ сделать это с помощью for loops и ifs:

import numpy as np
A = np.arange(1.,5.,1)
G = np.array([[1.,2.],[3.,4.]])

def calcF(G,A):
    N = A.size
    print A
    print N
    F = []
    for l in range(N):
        F.append(G/A[l])
        print F[l]
        for j in range(N):
            if j != l:
                F[l]*=((G - A[l])/(G + A[j]))*((A[l] - A[j])/(A[l] + A[j]))
    return F

F= calcF(G,A)
print F

Как для циклов, так и для операторов if относительно медленно, я ищу дурацкий остроумный способ сделать то же самое. У кого-нибудь есть идея?

1   2   2016-01-09 19:48:52

1 ответ:

Перечисленное в этом посте является векторизованным решением, делающим интенсивное использование мощной широковещательной функции NumPy после расширения размеров входных массивов до 3D и 4D случаев с np.newaxis/None в разных местах, согласно вычислениям, которые были задействованы. Вот реализация -

# Get size of A
N = A.size

# Perform "(G - A[l])/(G + A[j]))" in a vectorized manner
p1 = (G - A[:,None,None,None])/(G + A[:,None,None])

# Perform "((A[l] - A[j])/(A[l] + A[j]))" in a vectorized manner
p2 = ((A[:,None] - A)/(A[:,None] + A))

# Elementwise multiplications between the previously calculated parts
p3 = p1*p2[...,None,None]

# Set the escaped portion "j != l" output as "G/A[l]"
p3[np.eye(N,dtype=bool)] = G/A[:,None,None]
Fout = p3.prod(1)

# If you need separate arrays just like in the question, split it
Fout_split = np.array_split(Fout,N)

Пробный прогон -

In [284]: # Original inputs
     ...: A = np.arange(1.,5.,1)
     ...: G = np.array([[1.,2.],[3.,4.]])
     ...: 

In [285]: calcF(G,A)
Out[285]: 
[array([[-0.        , -0.00166667],
        [-0.01142857, -0.03214286]]), array([[-0.00027778,  0.        ],
        [ 0.00019841,  0.00126984]]), array([[  1.26984127e-03,   1.32275132e-04],
        [ -0.00000000e+00,  -7.93650794e-05]]), array([[-0.00803571, -0.00190476],
        [-0.00017857,  0.        ]])]

In [286]: vectorized_calcF(G,A) # Posted solution 
Out[286]: 
[array([[[-0.        , -0.00166667],
         [-0.01142857, -0.03214286]]]), array([[[-0.00027778,  0.        ],
         [ 0.00019841,  0.00126984]]]), array([[[  1.26984127e-03,   1.32275132e-04],
         [ -0.00000000e+00,  -7.93650794e-05]]]), array([[[-0.00803571, -0.00190476],
         [-0.00017857,  0.        ]]])]

Тест времени выполнения -

In [289]: # Larger inputs
     ...: A = np.random.randint(1,500,(400))
     ...: G = np.random.randint(1,400,(20,20))
     ...: 

In [290]: %timeit calcF(G,A)
1 loops, best of 3: 4.46 s per loop

In [291]: %timeit vectorized_calcF(G,A)  # Posted solution 
1 loops, best of 3: 1.87 s per loop

Векторизация с помощью NumPy / MATLAB: общий подход

Я чувствовал, что могу бросить свои два цента на мой общий подход, и я думаю, что другие следуют аналогичным стратегиям при попытке векторизации кодов, особенно в высокоуровневой платформе, такой как NumPy или MATLAB. Итак, вот краткий контрольный список вещей, которые можно было бы рассмотреть для Vectorization -

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

С чего начать векторизацию?Начните с самого глубокого (того цикла, где код повторяется больше всего) этапа вычисления и посмотрите, как можно расширить входные данные и ввести соответствующие вычисления. Внимательно следите за отслеживанием задействованных итераторов и соответствующим образом расширяйте измерения. Двигайтесь наружу на внешние петли, пока не будете удовлетворены выполненной векторизацией.

Как позаботиться об условных выражениях? для простых случаев, грубая сила вычисляет все и видит, как части IF / ELSE можно было бы позаботиться позже. Это было бы очень специфично для контекста.

Есть ли зависимости? если да, то посмотрите, можно ли проследить и реализовать зависимости соответствующим образом. Это может стать еще одной темой для обсуждения, но вот few examples я ввязался в это дело.