Как создать действительно вызываемый массив или матрицу в Python

Question

Я хотел бы взять матрицу, которая имеет все свои записи как функцию некоторой переменной x. Следовательно, B(x) будет выдавать N x N в идеале быстро. На самом деле это простая задача, если вы хотите напечатать матрицу с функциями в качестве записей. В качестве примера:

f1 = lambda x: 2*x
f2 = lambda x: x**2
B = lambda x : np.array([[f1(x),f2(x)],
                         [f2(x),f1(x)]])

Это наивно, так как не масштабируется в сценарии, где ваш массив большой и имеет множество функций. Распечатка этого заняла бы много времени для единственной проблемы. Обычно нужно создать пустой массив и использовать два цикла Python for, чтобы вычислить указанную функцию c для данной записи, а затем поместить выходные данные в массив. Затем массив возвращается.

Проблема с описанным выше методом заключается в том, что каждый раз, когда вызывается функция, она запускает их для циклов. Это замедляет работу, если вы хотите запустить функцию для набора данных со значениями x. Я пытался создать массив с возможностью вызова c, используя функцию lambdfiy Sympy. Для оценки x это кажется быстрее, чем для решения для l oop в чистом Python. Тем не менее, это ужасно перевешивает стоимость установки. Пожалуйста, смотрите мой код ниже для деталей.

Есть ли способ использовать функцию vectorize в Numpy, чтобы ускорить процесс? Сможете ли вы найти решение, которое быстрее, чем версия для l oop?

Я также играю с идеей (или называю это сном), где можно оценить весь набор данных X вместо каждого x отдельно. Как трансляция в Numpy. Например,

# Naive
result1 = [np.sin(x) for x in X]
# vs 
results2 = np.sin(X)

В любом случае, это довольно надумано. Вот код, который я написал. Пожалуйста, поиграйте с размером N, чтобы увидеть, насколько увлекательным является снижение скорости. Просто чтобы прояснить, я оценил свою программу в целом, и эта проблема вызываемого массива - то, где узкое место лежит.

import numpy as np
from sympy import symbols,lambdify,zeros
from time import time


def get_sympy(N):
    '''
    Creates a callable array using Sympys lambdfiy capabilites.
    This is truly a callable array.
    '''
    x = symbols('x')
    output = zeros(N,N)
    for i in range(N):
        for j in range(N):
            if i == j:
                output[i,j] = x**2
            elif i == 1:
                output[i,j] = x**3
            elif j == 0:
                output[i,j] = x**4
            else:
                output[i,j] = x
    return lambdify(x,output,'numpy')

def get_python(x,N):
    '''
    Uses Python loops to setup an array that mimic that of a callable array.
    It is not truly a callable array as the loops run on each call of 
    this function.
    '''
    output = np.zeros((N,N))
    f1 = lambda x: x**2
    f2 = lambda x: x**3
    f3 = lambda x: x**4
    for i in range(N):
        for j in range(N):
            if i == j:
                output[i,j] = f1(x)
            elif i == 1:
                output[i,j] = f2(x)
            elif j == 0:
                output[i,j] = f3(x)
            else:
                output[i,j] = x
    return output


if __name__ == '__main__':
    N = 30
    X = np.random.uniform()
    callable_sympy_array = get_sympy(N)
    callable_python_array = lambda x: get_python(x,N)
    t1 = time()
    sympy_result = callable_sympy_array(X)
    t2 = time()
    python_result = callable_python_array(X)
    t3 = time()
    sympy_func = get_sympy(N)
    t4 = time()
    sympy_time = t2-t1
    python_time = t3-t2
    sympy_setup_time = t4-t3

    print('\nSingle Evaluation Results:\n')
    print('Sympy: ',round(sympy_time,5))
    print('Python: ',round(python_time,5))
    print('Sympy + Setup',round(sympy_setup_time,5))

    evals = 100
    print('\nResults for {0} evaluations of a {1} by {1} array:\n'.format(evals,N))
    print('Sympy: ',sympy_setup_time + evals*sympy_time)
    print('Python: ',python_time*evals)

Answer 1

Быстрая numpy оценка требует применения встроенных скомпилированных операторов / функций для целых массивов. Любая итерация уровня python замедляет вас, как и оценка (общих) функций Python в скалярах. Быстрые вещи в основном ограничены операторами (например, **) и ufunc (np.sin, et c).

Ваша сгенерированная sympy функция иллюстрирует это:

В сеансе isympy:

In [65]: M = get_sympy(3)                                                                                 

с использованием самоанализа кода ipython:

In [66]: M??                                                                                              
Signature: M(x)
Docstring:
Created with lambdify. Signature:

func(x)

Expression:

Matrix([[x**2, x, x], [x**3, x**2, x**3], [x**4, x, x**2]])

Source code:

def _lambdifygenerated(x):
    return (array([[x**2, x, x], [x**3, x**2, x**3], [x**4, x, x**2]]))


Imported modules:
Source:   
def _lambdifygenerated(x):
    return (array([[x**2, x, x], [x**3, x**2, x**3], [x**4, x, x**2]]))
File:      /<lambdifygenerated-8>
Type:      function

Так что это функция в x с использованием операций numpy, ** оператор и создание массива. Точно так же, как если бы вы его набрали. sympy создает это с помощью лексических подстановок в своем символьном c коде, так что вы можете сказать, что он «печатает».

Он может работать на скаляре

In [67]: M(3)                                                                                             
Out[67]: 
array([[ 9,  3,  3],
       [27,  9, 27],
       [81,  3,  9]])

в массиве, здесь получается результат (3,3,3):

In [68]: M(np.arange(1,4))                                                                                
Out[68]: 
array([[[ 1,  4,  9],
        [ 1,  2,  3],
        [ 1,  2,  3]],

       [[ 1,  8, 27],
        [ 1,  4,  9],
        [ 1,  8, 27]],

       [[ 1, 16, 81],
        [ 1,  2,  3],
        [ 1,  4,  9]]])

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

Ваш get_python не будет принимать массив x, прежде всего потому, что

 output = np.zeros((N,N))

имеет фиксированный размер; использование np.zeros((N,N)+x.shape), x.dtype) может обойти это.

В любом случае, оно будет медленным из-за итерации уровня python при каждом вызове.

===

Было бы быстрее, если бы вы попытались назначить группы элементов. Например, в этом случае:

In [76]: output = np.zeros((3,3),int)                                                                     
In [77]: output[:] = 3                                                                                    
In [78]: output[:,0]=3**4                                                                                 
In [79]: output[1,:]=3**3                                                                                 
In [80]: output[np.arange(3),np.arange(3)]=3**2                                                           

In [81]: output                                                                                           
Out[81]: 
array([[ 9,  3,  3],
       [27,  9, 27],
       [81,  3,  9]])

===

frompyfunc - удобный инструмент для подобных случаев. В некоторых случаях он предлагает увеличение скорости в 2 раза по сравнению с прямой итерацией. Но даже без этого он может сделать код более кратким.

Например, быстрое описание вашего примера:

In [82]: def foo(x,i,j): 
    ...:     if i==j: return x**2 
    ...:     if i==1: return x**3 
    ...:     if j==0: return x**4 
    ...:     return x                                                                                             
In [83]: f = np.frompyfunc(foo, 3, 1)                                                                     

In [84]: f(3,np.arange(3)[:,None], np.arange(3))                                                          
Out[84]: 
array([[9, 3, 3],
       [27, 9, 27],
       [81, 3, 9]], dtype=object)

и для примера Out[68]:

In [98]: f(np.arange(1,4),np.arange(3)[:,None,None], np.arange(3)[:,None]).shape                          
Out[98]: (3, 3, 3)

In [99]: timeit f(np.arange(1,4),np.arange(3)[:,None,None], np.arange(3)[:,None]).shape                   
23 µs ± 471 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [100]: timeit M(np.arange(1,4))                                                                        
21.7 µs ± 440 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

Оценивается по скаляру x, мой f примерно такой же скорости, как ваш get_python.

In [115]: MM = get_sympy(30)                                                                              
In [116]: timeit MM(3)                                                                                    
109 µs ± 112 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [117]: timeit get_python(3,30)                                                                         
241 µs ± 2.06 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

In [118]: timeit f(3,np.arange(30)[:,None], np.arange(30)).astype(int)                                    
254 µs ± 1.2 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

Answer 2

Мне очень нравится принятый ответ. Я просто хотел бы поделиться своим собственным решением.

Примером вызываемых матриц, которые я использовал, были игрушечные примеры. Фактически вызываемые матрицы - это цепи Маркова, для которых к их строкам применено преобразование Lo git. Я также буду брать производные этих матриц по некоторым параметрам. В любом случае, дело в том, что я делаю все это через Sympy, а не вручную. Следовательно, имеет смысл использовать функцию lambdify для моих результатов. Дополнительным преимуществом является то, что получается довольно быстрая вызываемая матрица.

Я просто хотел бы отметить, что вычисления и процесс lambdify требуют больших вычислительных ресурсов. Но у нас есть cloudpickle как наш друг.

Итак, что я сделал, так это вычислил свои вызываемые матрицы для каждого измерения N в range(2,500) и затем сохранил их в большом словаре, который сериализован. Примечание: cloudpickle очень устойчив и является единственным из pickle, cpickle или dill, чтобы справиться с этим без каких-либо дополнительных требований к настройке.

Вот небольшой пример:

from cloudpickle import dump,load

callable_arrays = dict()
for N in range(2,500):
    callable_arrays[N] = get_sympy(N)

# Serialize the dictionary
with open('callable_array_file','wb') as file:
    dump(callable_arrays,file)

# We can write a re-usable function to access callable arrays
def get_callable_array(N):
    output = None
    with open('callable_array_file','rb') as file:
        output = load(file)[N]
    return output

Возможно, это можно уточнить, но я доволен этой идеей. Сюрпризом дня стало то, что Sympy может генерировать вызываемый массив. Принятый ответ, данный @hpaulji, подробно показывает, почему это так.