Poprawa wydajności w transformacji tablic numpy

Oct 31 2020

Biorąc pod uwagę trzy numpytablice 1D, chcę je przekształcić w następujący sposób:

import numpy as np

Xd = np.asarray([0, 0,   1,   1,   0.5])
Yd = np.asarray([0, 0,   0,   2.5, 2.5])
Zd = np.asarray([0, 1.5, 1.5, 1.5, 1.5])

points = np.stack([Xd, Yd, Zd], axis=1).reshape(-1, 1, 3)
segments = np.concatenate([points[:-1], points[1:]], axis = 1)    

print(segments.shape)
print(segments)

Wynik:

(4, 2, 3)
[[[0.  0.  0. ]
  [0.  0.  1.5]]

 [[0.  0.  1.5]
  [1.  0.  1.5]]

 [[1.  0.  1.5]
  [1.  2.5 1.5]]

 [[1.  2.5 1.5]
  [0.5 2.5 1.5]]]

Czy istnieje sposób na poprawę wydajności tej transformacji?

tło

Ta transformacja jest konieczna, aby użyć XYZwspółrzędnych w matplotlibz Line3DCollection. Do tej pory widziałem tylko odmiany powyższego kodu, ale przy tysiącach współrzędnych lub danych interpolowanych dla lepszej rozdzielczości konieczne jest zoptymalizowane podejście.

Podsumowanie

Dzięki @Mercury można stwierdzić, że w przypadku krótszych tablic (<1k długości) odpowiedź @Miguel działa lepiej, ale podejście @mathfux skaluje się znacznie lepiej, gdy tablice stają się dłuższe.

Odpowiedzi

3 mathfux Oct 31 2020 at 13:03

Wygląda na to, że próbujesz obrócić okno kształtu (2, 3)w tablicy 2D. Jest to podobne do splotu obrazu, który można wykonać np.lib.stride_tricksw bardzo wydajny sposób.

a = np.transpose([Xd, Yd, Zd])
window = (2, 3)
view_shape = (len(a) - window[0] + 1,) + window # (4,2,3) if len(a) == 5
sub_matrix = np.lib.stride_tricks.as_strided(a, shape = view_shape, strides = (a.itemsize,) + a.strides)
>>> sub_matrix
array([[[0. , 0. , 0. ],
        [0. , 0. , 1.5]],

       [[0. , 0. , 1.5],
        [1. , 0. , 1.5]],

       [[1. , 0. , 1.5],
        [1. , 2.5, 1.5]],

       [[1. , 2.5, 1.5],
        [0.5, 2.5, 1.5]]])

Zauważ, że np.lib.stride_tricksjest to bardzo wydajne w porównaniu z innymi sposobami.

4 Miguel Oct 31 2020 at 12:22

Ogólna rada, jeśli zależy Ci na szybkości, powinna generalnie unikać łączenia i łączenia, ponieważ często oznacza to wielokrotne kopiowanie tych samych danych.

W każdym razie, oto jak bym to zrobił, nieco dłuższy kod, ale nie wykonuje więcej pracy niż potrzeba

n = len(Xd)
segments = np.empty((n-1, 2, 3))

segments[:,0,0] = Xd[:-1]
segments[:,1,0] = Xd[1:]

segments[:,0,1] = Yd[:-1]
segments[:,1,1] = Yd[1:]

segments[:,0,2] = Zd[:-1]
segments[:,1,2] = Zd[1:]

[EDIT] - Następujące powstał dla nauki / zabawy, nie reproduc e

Więc spróbowałem sprawdzić, czy mógłbym wycisnąć trochę więcej wydajności z odpowiedzi @mathfux i wyszedłem z tym brzydkim kodem:

a = np.empty(3*n)
a[:n]    = Xd
a[n:n+n] = Yd
a[n+n:]  = Zd

interface = dict(a.__array_interface__)
interface['shape'] = (n-1, 2, 3)
interface['strides'] = (a.itemsize, a.itemsize, n*a.itemsize)
segments= np.array(np.lib.stride_tricks.DummyArray(interface, base=a), copy=False)

Na moim komputerze jest to mierzalnie szybsze (do ~ 30% w zależności od wielkości wkładu). Zyski są częściowo spowodowane budową ai pomijaniem kontrolias_strided

1 Mercury Oct 31 2020 at 13:03

Oto kilka testów czasowych na większych tablicach, co wyjaśnia różnicę.

import numpy as np
from timeit import timeit

# original
def f1(x, y, z):
    points = np.stack([x, y, z], axis=1).reshape(-1, 1, 3)
    return np.concatenate([points[:-1], points[1:]], axis = 1)

# preallocating and then assigning
def f2(x, y, z):
    segments = np.empty((len(x)-1, 2, 3))

    segments[:,0,0] = x[:-1]
    segments[:,1,0] = x[1:]

    segments[:,0,1] = y[:-1]
    segments[:,1,1] = y[1:]

    segments[:,0,2] = z[:-1]
    segments[:,1,2] = z[1:]
    return segments

# stacking, but in one go
def f3(x, y, z):
    segments = np.stack([x[:-1], y[:-1], z[:-1], x[1:], y[1:],z[1:]], axis=1)
    return segments.reshape(-1, 2, 3)

# list comparison
def f4(x, y, z):
    z_ = [i for i in zip(x,y,z)]
    return [[[z_[i]],[z_[i+1]]] for i in range(len(z_)-1)]

#np.lib.stride_tricks approach
def f5(x, y, z):
    a = np.transpose([x, y, z])
    window = (2, 3)
    view_shape = (len(a) - window[0] + 1,) + window # (4,2,3) if len(a) == 5
    return np.lib.stride_tricks.as_strided(a, shape = view_shape, strides = (a.itemsize,) + a.strides)
    

ntime = 5000 #number of test runs
nxd = 500    #array length

Xd = np.random.randn(nxd)
Yd = np.random.randn(nxd)
Zd = np.random.randn(nxd)

print(timeit(lambda: f1(Xd, Yd, Zd), number=ntime))
#0.11369249999999999

print(timeit(lambda: f2(Xd, Yd, Zd), number=ntime))
#0.0480651

print(timeit(lambda: f3(Xd, Yd, Zd), number=ntime))
#0.10202380000000003

print(timeit(lambda: f4(Xd, Yd, Zd), number=ntime))
#1.8407391

print(timeit(lambda: f5(Xd, Yd, Zd), number=ntime))
#0.09132560000000023
    
ntime = 50     #number of test runs
nxd = 500000   #array length

Xd = np.random.randn(nxd)
Yd = np.random.randn(nxd)
Zd = np.random.randn(nxd)

print(timeit(lambda: f1(Xd, Yd, Zd), number=ntime))
#1.7519548999999999

print(timeit(lambda: f2(Xd, Yd, Zd), number=ntime))
#1.504727

print(timeit(lambda: f3(Xd, Yd, Zd), number=ntime))
#1.5010566

print(timeit(lambda: f4(Xd, Yd, Zd), number=ntime))
#22.6208157

print(timeit(lambda: f5(Xd, Yd, Zd), number=ntime))
#0.46465339999999955

Jak widać, droga @ Miguela jest drogą: wstępne przydzielenie tablicy, a następnie przypisanie jest najbardziej wydajnym sposobem na zrobienie tego. Nawet jeśli ułożysz je w sprytniejszy sposób, jak w f3 (), nadal jest wolniejszy niż f2 (). Ale nic nie przebije f5 (), gdy długość tablicy znacznie się zwiększy.

venky__ Oct 31 2020 at 12:31

Okazało się, że jest to szybsze niż kod @ Miguel.

z = [i for i in zip(Xd,Yd,Zd)]
segments = [[[z[i]],[z[i+1]]] for i in range(len(z)-1)]