Używanie numpy do zbudowania tablicy wszystkich kombinacji dwóch tablic

143

Próbuję przeanalizować przestrzeń parametrów funkcji 6-parametrowej, aby zbadać jej zachowanie numeryczne, zanim spróbuję zrobić z nią coś złożonego, więc szukam skutecznego sposobu, aby to zrobić.

Moja funkcja przyjmuje wartości zmiennoprzecinkowe, które jako dane wejściowe przyjmują 6-dim tablicę numpy. Na początku próbowałem zrobić tak:

Najpierw stworzyłem funkcję, która pobiera 2 tablice i generuje tablicę ze wszystkimi kombinacjami wartości z dwóch tablic

from numpy import *
def comb(a,b):
    c = []
    for i in a:
        for j in b:
            c.append(r_[i,j])
    return c

Następnie reduce()stosowałem to do m kopii tej samej tablicy:

def combs(a,m):
    return reduce(comb,[a]*m)

A potem oceniam moją funkcję w ten sposób:

values = combs(np.arange(0,1,0.1),6)
for val in values:
    print F(val)

To działa, ale jest zbyt wolne. Wiem, że przestrzeń parametrów jest ogromna, ale nie powinno to być takie wolne. W tym przykładzie pobrałem tylko 10 6 (milion) punktów, a samo utworzenie tablicy zajęło więcej niż 15 sekund values.

Czy znasz lepszy sposób na zrobienie tego z numpy?

Mogę zmodyfikować sposób, w jaki funkcja Fprzyjmuje argumenty, jeśli jest to konieczne.

Rafael S. Calsaverini
źródło
Aby znaleźć najszybszy produkt kartezjański, jaki znalazłem, zobacz tę odpowiedź . (Ponieważ pytanie jest sformułowane zupełnie inaczej niż to, uważam, że pytania nie są duplikatami, ale najlepszym rozwiązaniem jest to samo).
senderle

Odpowiedzi:

127

W nowszej wersji numpy(> 1.8.x) numpy.meshgrid()zapewnia znacznie szybszą implementację:

Rozwiązanie @ pv

In [113]:

%timeit cartesian(([1, 2, 3], [4, 5], [6, 7]))
10000 loops, best of 3: 135 µs per loop
In [114]:

cartesian(([1, 2, 3], [4, 5], [6, 7]))

Out[114]:
array([[1, 4, 6],
       [1, 4, 7],
       [1, 5, 6],
       [1, 5, 7],
       [2, 4, 6],
       [2, 4, 7],
       [2, 5, 6],
       [2, 5, 7],
       [3, 4, 6],
       [3, 4, 7],
       [3, 5, 6],
       [3, 5, 7]])

numpy.meshgrid()Kiedyś był tylko 2D, teraz jest zdolny do ND. W tym przypadku 3D:

In [115]:

%timeit np.array(np.meshgrid([1, 2, 3], [4, 5], [6, 7])).T.reshape(-1,3)
10000 loops, best of 3: 74.1 µs per loop
In [116]:

np.array(np.meshgrid([1, 2, 3], [4, 5], [6, 7])).T.reshape(-1,3)

Out[116]:
array([[1, 4, 6],
       [1, 5, 6],
       [2, 4, 6],
       [2, 5, 6],
       [3, 4, 6],
       [3, 5, 6],
       [1, 4, 7],
       [1, 5, 7],
       [2, 4, 7],
       [2, 5, 7],
       [3, 4, 7],
       [3, 5, 7]])

Zwróć uwagę, że kolejność wyniku końcowego jest nieco inna.

CT Zhu
źródło
14
np.stack(np.meshgrid([1, 2, 3], [4, 5], [6, 7]), -1).reshape(-1, 3)wyda właściwą kolejność
Eric
@CT Zhu Czy istnieje łatwy sposób na przekształcenie tego, aby zamiast tego jako dane wejściowe używana była macierz przechowująca różne tablice jako kolumny?
Dole
2
Należy zauważyć, że meshgrid działa tylko dla mniejszych zestawów zakresów, mam duży i otrzymuję błąd: ValueError: maksymalny obsługiwany wymiar ndarray to 32, znaleziono 69
mikkom
158

Oto prosta implementacja. To około 5 razy szybciej niż przy użyciu itertools.


import numpy as np

def cartesian(arrays, out=None):
    """
    Generate a cartesian product of input arrays.

    Parameters
    ----------
    arrays : list of array-like
        1-D arrays to form the cartesian product of.
    out : ndarray
        Array to place the cartesian product in.

    Returns
    -------
    out : ndarray
        2-D array of shape (M, len(arrays)) containing cartesian products
        formed of input arrays.

    Examples
    --------
    >>> cartesian(([1, 2, 3], [4, 5], [6, 7]))
    array([[1, 4, 6],
           [1, 4, 7],
           [1, 5, 6],
           [1, 5, 7],
           [2, 4, 6],
           [2, 4, 7],
           [2, 5, 6],
           [2, 5, 7],
           [3, 4, 6],
           [3, 4, 7],
           [3, 5, 6],
           [3, 5, 7]])

    """

    arrays = [np.asarray(x) for x in arrays]
    dtype = arrays[0].dtype

    n = np.prod([x.size for x in arrays])
    if out is None:
        out = np.zeros([n, len(arrays)], dtype=dtype)

    m = n / arrays[0].size
    out[:,0] = np.repeat(arrays[0], m)
    if arrays[1:]:
        cartesian(arrays[1:], out=out[0:m, 1:])
        for j in xrange(1, arrays[0].size):
            out[j*m:(j+1)*m, 1:] = out[0:m, 1:]
    return out
pv.
źródło
47
kiedykolwiek rozważałeś przesłanie tego do numpy? to nie pierwszy raz, kiedy szukałem tej funkcji i znalazłem Twój post.
endolith
1
W tej implementacji jest błąd. Na przykład tablice łańcuchów: tablice [0] .dtype = "| S3" i tablice [1] .dtype = "| S5". Istnieje więc potrzeba znalezienia najdłuższego łańcucha wejściowego i użycia jego typu w out = np.zeros ([n, len (tablice)], dtype = dtype)
norecces
38
Do Twojej wiadomości: wydaje się, że znalazło się w pakiecie scikit-learn wfrom sklearn.utils.extmath import cartesian
Gus
2
Właśnie sobie uświadomiłem: to trochę różni się od itertools.combinations, ponieważ ta funkcja przestrzega kolejności wartości, podczas gdy kombinacje nie, więc ta funkcja zwraca więcej wartości niż kombinacje. Wciąż imponujące, ale niestety nie to, czego szukałem :(
David Marx
6
TypeError: slice indices must be integers or None or have an __index__ methodwyrzucony przezcartesian(arrays[1:], out=out[0:m,1:])
Boern
36

itertools.combinations to generalnie najszybszy sposób na uzyskanie kombinacji z kontenera Pythona (jeśli faktycznie potrzebujesz kombinacji, tj. układów BEZ powtórzeń i niezależnych od kolejności; wydaje się, że nie to robi twój kod, ale nie mogę powiedz, czy dzieje się tak dlatego, że twój kod zawiera błędy, czy też dlatego, że używasz złej terminologii).

Jeśli chcesz czegoś innego niż kombinacje, być może inne iteratory w itertools productlub permutationsmogą ci służyć lepiej. Na przykład wygląda na to, że Twój kod jest mniej więcej taki sam, jak:

for val in itertools.product(np.arange(0, 1, 0.1), repeat=6):
    print F(val)

Wszystkie te iteratory dają krotki, a nie listy lub tablice numpy, więc jeśli twój F jest wybredny, jeśli chodzi o uzyskanie konkretnej tablicy numpy, będziesz musiał zaakceptować dodatkowe obciążenie związane z konstruowaniem lub czyszczeniem i ponownym wypełnianiem jednej na każdym kroku.

Alex Martelli
źródło
8

Możesz zrobić coś takiego

import numpy as np

def cartesian_coord(*arrays):
    grid = np.meshgrid(*arrays)        
    coord_list = [entry.ravel() for entry in grid]
    points = np.vstack(coord_list).T
    return points

a = np.arange(4)  # fake data
print(cartesian_coord(*6*[a])

co daje

array([[0, 0, 0, 0, 0, 0],
   [0, 0, 0, 0, 0, 1],
   [0, 0, 0, 0, 0, 2],
   ..., 
   [3, 3, 3, 3, 3, 1],
   [3, 3, 3, 3, 3, 2],
   [3, 3, 3, 3, 3, 3]])
felippe
źródło
2
Czy istnieje sposób, aby NumPy akceptował więcej niż 32 tablice dla siatki mesh? Ta metoda działa dla mnie, o ile nie przekażę więcej niż 32 tablice.
Joelmob
8

Następująca implementacja numpy powinna wynosić ok. 2x szybkość udzielonej odpowiedzi:

def cartesian2(arrays):
    arrays = [np.asarray(a) for a in arrays]
    shape = (len(x) for x in arrays)

    ix = np.indices(shape, dtype=int)
    ix = ix.reshape(len(arrays), -1).T

    for n, arr in enumerate(arrays):
        ix[:, n] = arrays[n][ix[:, n]]

    return ix
Stefan van der Walt
źródło
1
Wygląda dobrze. Z moich podstawowych testów wynika, że ​​wygląda to szybciej niż pierwotna odpowiedź dla wszystkich par, trójek i 4 krotek {1, 2, ..., 100}. Potem wygrywa oryginalna odpowiedź. Również dla przyszłych czytelników, którzy chcą wygenerować wszystkie k-krotki {1, ..., n}, np.indices((n,...,n)).reshape(k,-1).Twystarczą.
jme
Działa to tylko dla liczb całkowitych, a zaakceptowana odpowiedź działa również dla liczb zmiennoprzecinkowych.
FJC
7

Wygląda na to, że potrzebujesz siatki do oceny funkcji, w takim przypadku możesz użyć numpy.ogrid(otwórz) lub numpy.mgrid(uzupełnij):

import numpy
my_grid = numpy.mgrid[[slice(0,1,0.1)]*6]
Steabert
źródło
6

możesz użyć np.array(itertools.product(a, b))

William Song
źródło
np.array (list (itertools.product (l, l2)))
ZirconCode
4

Oto jeszcze jeden sposób, używając czystego NumPy, bez rekursji, bez rozumienia list i bez jawnych pętli for. Jest około 20% wolniejsze niż oryginalna odpowiedź i opiera się na np.meshgrid.

def cartesian(*arrays):
    mesh = np.meshgrid(*arrays)  # standard numpy meshgrid
    dim = len(mesh)  # number of dimensions
    elements = mesh[0].size  # number of elements, any index will do
    flat = np.concatenate(mesh).ravel()  # flatten the whole meshgrid
    reshape = np.reshape(flat, (dim, elements)).T  # reshape and transpose
    return reshape

Na przykład,

x = np.arange(3)
a = cartesian(x, x, x, x, x)
print(a)

daje

[[0 0 0 0 0]
 [0 0 0 0 1]
 [0 0 0 0 2]
 ..., 
 [2 2 2 2 0]
 [2 2 2 2 1]
 [2 2 2 2 2]]
étale-cohomology
źródło
3

Aby uzyskać czystą numpy implementację iloczynu kartezjańskiego tablic 1D (lub płaskich list Pythona), po prostu użyj meshgrid(), obróć osie transpose()i zmień kształt na pożądany wynik:

 def cartprod(*arrays):
     N = len(arrays)
     return transpose(meshgrid(*arrays, indexing='ij'), 
                      roll(arange(N + 1), -1)).reshape(-1, N)

Zauważ, że ma to konwencję najszybszej zmiany ostatniej osi („styl C” lub „wiersz główny”).

In [88]: cartprod([1,2,3], [4,8], [100, 200, 300, 400], [-5, -4])
Out[88]: 
array([[  1,   4, 100,  -5],
       [  1,   4, 100,  -4],
       [  1,   4, 200,  -5],
       [  1,   4, 200,  -4],
       [  1,   4, 300,  -5],
       [  1,   4, 300,  -4],
       [  1,   4, 400,  -5],
       [  1,   4, 400,  -4],
       [  1,   8, 100,  -5],
       [  1,   8, 100,  -4],
       [  1,   8, 200,  -5],
       [  1,   8, 200,  -4],
       [  1,   8, 300,  -5],
       [  1,   8, 300,  -4],
       [  1,   8, 400,  -5],
       [  1,   8, 400,  -4],
       [  2,   4, 100,  -5],
       [  2,   4, 100,  -4],
       [  2,   4, 200,  -5],
       [  2,   4, 200,  -4],
       [  2,   4, 300,  -5],
       [  2,   4, 300,  -4],
       [  2,   4, 400,  -5],
       [  2,   4, 400,  -4],
       [  2,   8, 100,  -5],
       [  2,   8, 100,  -4],
       [  2,   8, 200,  -5],
       [  2,   8, 200,  -4],
       [  2,   8, 300,  -5],
       [  2,   8, 300,  -4],
       [  2,   8, 400,  -5],
       [  2,   8, 400,  -4],
       [  3,   4, 100,  -5],
       [  3,   4, 100,  -4],
       [  3,   4, 200,  -5],
       [  3,   4, 200,  -4],
       [  3,   4, 300,  -5],
       [  3,   4, 300,  -4],
       [  3,   4, 400,  -5],
       [  3,   4, 400,  -4],
       [  3,   8, 100,  -5],
       [  3,   8, 100,  -4],
       [  3,   8, 200,  -5],
       [  3,   8, 200,  -4],
       [  3,   8, 300,  -5],
       [  3,   8, 300,  -4],
       [  3,   8, 400,  -5],
       [  3,   8, 400,  -4]])

Jeśli chcesz zmienić pierwszą oś najszybciej („styl FORTRAN” lub „kolumna główna”), po prostu zmień orderparametr w reshape()ten sposób:reshape((-1, N), order='F')

RBF06
źródło
1

Pandy mergeoferują naiwne, szybkie rozwiązanie problemu:

# given the lists
x, y, z = [1, 2, 3], [4, 5], [6, 7]

# get dfs with same, constant index 
x = pd.DataFrame({'x': x}, index=np.repeat(0, len(x))
y = pd.DataFrame({'y': y}, index=np.repeat(0, len(y))
z = pd.DataFrame({'z': z}, index=np.repeat(0, len(z))

# get all permutations stored in a new df
df = pd.merge(x, pd.merge(y, z, left_index=True, righ_index=True),
              left_index=True, right_index=True)
simone
źródło