Evaluasi fungsi yang efisien di setiap sel dari array NumPy

Diberikan array NumPy A , apa cara tercepat / paling efisien untuk menerapkan fungsi yang sama , f , ke setiap sel?

1. Misalkan kita akan tetapkan untuk A (i, j) yang f (A (i, j)) .

2. Fungsi, f , tidak memiliki keluaran biner, sehingga operasi mask tidak akan membantu.

Apakah iterasi loop ganda yang "jelas" (melalui setiap sel) merupakan solusi yang optimal?

Anda bisa melakukan vektorisasi fungsi dan kemudian menerapkannya langsung ke array Numpy setiap kali Anda membutuhkannya:

import numpy as np

def f(x):
    return x * x + 3 * x - 2 if x > 0 else x * 5 + 8

f = np.vectorize(f)  # or use a different name if you want to keep the original f

result_array = f(A)  # if A is your Numpy array

Mungkin lebih baik untuk menentukan tipe keluaran eksplisit secara langsung saat melakukan vektorisasi:

f = np.vectorize(f, otypes=[np.float])

Pertanyaan serupa adalah: Memetakan array NumPy di ​​tempat . Jika Anda dapat menemukan ufunc untuk f () Anda, maka Anda harus menggunakan parameter out.

Jika Anda bekerja dengan angka dan f(A(i,j)) = f(A(j,i)), Anda bisa menggunakan scipy.spatial.distance.cdist yang mendefinisikan f sebagai jarak antara A(i)dan A(j).

Saya yakin saya telah menemukan solusi yang lebih baik. Ide untuk mengubah fungsi menjadi fungsi universal python (lihat dokumentasi ), yang dapat melakukan komputasi paralel di bawah tenda.

Seseorang dapat menulis sendiri disesuaikan ufuncdalam C, yang pasti lebih efisien, atau dengan memanggil np.frompyfunc, yang merupakan metode pabrik bawaan. Setelah pengujian, ini lebih efisien daripada np.vectorize:

f = lambda x, y: x * y
f_arr = np.frompyfunc(f, 2, 1)
vf = np.vectorize(f)
arr = np.linspace(0, 1, 10000)

%timeit f_arr(arr, arr) # 307ms
%timeit f_arr(arr, arr) # 450ms

Saya juga telah menguji sampel yang lebih besar, dan peningkatannya proporsional. Untuk perbandingan kinerja metode lain, lihat posting ini

Ketika 2d-array (atau nd-array) adalah C- atau F-contiguous, maka tugas memetakan sebuah fungsi ke dalam array 2d secara praktis sama dengan tugas memetakan fungsi ke dalam array 1d - kita hanya harus melihatnya seperti itu, mis np.ravel(A,'K'). via .

Solusi yang memungkinkan untuk 1d-array telah dibahas misalnya di sini .

Namun, ketika memori dari 2d-array tidak bersebelahan, maka situasinya sedikit lebih rumit, karena seseorang ingin menghindari kemungkinan cache miss jika sumbu ditangani dalam urutan yang salah.

Numpy sudah memiliki mesin untuk memproses sumbu dalam urutan terbaik. Salah satu kemungkinan untuk menggunakan mesin ini adalah np.vectorize. Namun, dokumentasi numpy np.vectorizemenyatakan bahwa itu "disediakan terutama untuk kenyamanan, bukan untuk kinerja" - fungsi python yang lambat tetap menjadi fungsi python yang lambat dengan seluruh overhead terkait! Masalah lainnya adalah konsumsi memori yang besar - lihat contoh SO-post ini .

Ketika seseorang ingin memiliki kinerja fungsi-C tetapi menggunakan mesin numpy, solusi yang baik adalah menggunakan numba untuk pembuatan ufunc, misalnya:

# runtime generated C-function as ufunc
import numba as nb
@nb.vectorize(target="cpu")
def nb_vf(x):
    return x+2*x*x+4*x*x*x

Ini mudah berdetak np.vectorizetetapi juga ketika fungsi yang sama akan dilakukan sebagai perkalian / penambahan numpy-array, yaitu

# numpy-functionality
def f(x):
    return x+2*x*x+4*x*x*x

# python-function as ufunc
import numpy as np
vf=np.vectorize(f)
vf.__name__="vf"

Lihat lampiran jawaban ini untuk kode pengukuran waktu:

Versi Numba (hijau) sekitar 100 kali lebih cepat daripada fungsi python (yaitu np.vectorize), yang tidak mengherankan. Tetapi juga sekitar 10 kali lebih cepat daripada fungsionalitas numpy, karena versi numbas tidak memerlukan array perantara dan karenanya menggunakan cache lebih efisien.

Meskipun pendekatan ufunc numba adalah pertukaran yang baik antara kegunaan dan kinerja, itu masih bukan yang terbaik yang bisa kami lakukan. Namun tidak ada solusi terbaik atau pendekatan terbaik untuk tugas apa pun - orang harus memahami apa saja batasannya dan bagaimana hal itu dapat dikurangi.

Misalnya, untuk fungsi transendental (misalnya exp, sin, cos) Numba tidak memberikan keuntungan apa pun atas numpy ini np.exp(tidak ada array sementara yang dibuat - sumber utama kecepatan-up). Namun, instalasi Anaconda saya menggunakan VML Intel untuk vektor yang lebih besar dari 8192 - tidak dapat melakukannya jika memori tidak berdekatan. Jadi mungkin lebih baik untuk menyalin elemen ke memori yang berdekatan agar dapat menggunakan VML Intel:

import numba as nb
@nb.vectorize(target="cpu")
def nb_vexp(x):
    return np.exp(x)

def np_copy_exp(x):
    copy = np.ravel(x, 'K')
    return np.exp(copy).reshape(x.shape) 

Untuk keadilan perbandingan, saya telah mematikan paralelisasi VML (lihat kode di lampiran):

Seperti yang bisa dilihat, begitu VML dijalankan, overhead penyalinan lebih dari kompensasi. Namun begitu data menjadi terlalu besar untuk cache L3, keuntungannya minimal karena tugas sekali lagi terikat pada memori-bandwidth.

Di sisi lain, numba juga dapat menggunakan SVML Intel, seperti yang dijelaskan dalam posting ini :

from llvmlite import binding
# set before import
binding.set_option('SVML', '-vector-library=SVML')

import numba as nb

@nb.vectorize(target="cpu")
def nb_vexp_svml(x):
    return np.exp(x)

dan menggunakan VML dengan hasil paralelisasi:

Versi numba memiliki overhead yang lebih sedikit, tetapi untuk beberapa ukuran VML mengalahkan SVML meskipun ada overhead penyalinan tambahan - yang tidak mengejutkan karena ufunc numba tidak diparalelkan.

Daftar:

A. perbandingan fungsi polinom:

import perfplot
perfplot.show(
    setup=lambda n: np.random.rand(n,n)[::2,::2],
    n_range=[2**k for k in range(0,12)],
    kernels=[
        f,
        vf, 
        nb_vf
        ],
    logx=True,
    logy=True,
    xlabel='len(x)'
    ) 

B. perbandingan dari exp:

import perfplot
import numexpr as ne # using ne is the easiest way to set vml_num_threads
ne.set_vml_num_threads(1)
perfplot.show(
    setup=lambda n: np.random.rand(n,n)[::2,::2],
    n_range=[2**k for k in range(0,12)],
    kernels=[
        nb_vexp, 
        np.exp,
        np_copy_exp,
        ],
    logx=True,
    logy=True,
    xlabel='len(x)',
    )

Semua jawaban di atas sebanding dengan baik, tetapi jika Anda perlu menggunakan fungsi kustom untuk pemetaan, dan Anda punya numpy.ndarray, dan Anda perlu mempertahankan bentuk array.

Saya telah membandingkan hanya dua, tetapi itu akan mempertahankan bentuk ndarray. Saya telah menggunakan array dengan 1 juta entri untuk perbandingan. Di sini saya menggunakan fungsi persegi. Saya menyajikan kasus umum untuk array dimensi n. Untuk dua dimensi buat saja iter2D.

import numpy, time

def A(e):
    return e * e

def timeit():
    y = numpy.arange(1000000)
    now = time.time()
    numpy.array([A(x) for x in y.reshape(-1)]).reshape(y.shape)        
    print(time.time() - now)
    now = time.time()
    numpy.fromiter((A(x) for x in y.reshape(-1)), y.dtype).reshape(y.shape)
    print(time.time() - now)
    now = time.time()
    numpy.square(y)  
    print(time.time() - now)

Keluaran

>>> timeit()
1.162431240081787    # list comprehension and then building numpy array
1.0775556564331055   # from numpy.fromiter
0.002948284149169922 # using inbuilt function

di sini Anda dapat dengan jelas melihat numpy.fromiterfungsi persegi pengguna, gunakan pilihan Anda. Jika fungsi Anda bergantung pada i, j itu adalah indeks array, iterasi pada ukuran array seperti for ind in range(arr.size), gunakan numpy.unravel_indexuntuk mendapatkan i, j, ..berdasarkan indeks 1D Anda dan bentuk array numpy.unravel_index

Jawaban ini terinspirasi oleh jawaban saya atas pertanyaan lain di sini