Utilisation de numpy pour construire un tableau de toutes les combinaisons de deux tableaux
j'essaie de parcourir l'espace paramètres d'une fonction à 6 paramètres pour étudier son comportement numérique avant d'essayer de faire quoi que ce soit de complexe avec elle donc je cherche un moyen efficace de le faire.
ma fonction prend en entrée des valeurs float données par un tableau de numpy 6-dim. Ce que j'ai essayé de faire au départ était ceci:
J'ai d'abord créé une fonction qui prend 2 Tableaux et génère un tableau avec toutes les combinaisons de valeurs des deux tableaux
from numpy import *
def comb(a,b):
c = []
for i in a:
for j in b:
c.append(r_[i,j])
return c
puis j'ai utilisé reduce() pour appliquer cela à m copies du même tableau:
def combs(a,m):
return reduce(comb,[a]*m)
et puis j'évalue ma fonction comme ceci:
values = combs(np.arange(0,1,0.1),6)
for val in values:
print F(val)
ça marche mais c'est trop lent. Je sais que l'espace des paramètres est énorme, mais cela ne devrait pas être si lent. J'ai seulement échantillonné 10 6 (un million) points dans cet exemple et il a fallu plus de 15 secondes juste pour créer le tableau values .
connais-tu un moyen plus efficace de le faire avec numpy?
je peux modifier la façon la fonction F prend ses arguments si c'est nécessaire.
9 réponses
dans la nouvelle version de numpy (>1.8.x), np.meshgrid fournit une mise en œuvre beaucoup plus rapide:
@pv de la solution
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 utiliser pour être 2D seulement, maintenant il est capable de ND. Dans ce cas, 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]])
Notez que l'ordre de la résultante finale est légèrement différente.
Voici une implémentation pure et simple. C'est ca. 5× plus rapide que l'utilisation d'itértools.
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
itertools.combinaisons est en général le moyen le plus rapide pour obtenir des combinaisons à partir d'un conteneur Python (Si vous voulez en fait des combinaisons, c'est-à-dire des arrangements sans répétitions et indépendamment de l'ordre; ce n'est pas ce que votre code semble faire, mais je ne peux pas dire si c'est parce que votre code est bogué ou parce que vous utilisez la mauvaise terminologie).
si vous voulez quelque chose de différent que des combinaisons peut-être d'autres itérateurs dans itertools, product ou permutations , pourrait mieux vous servir. Par exemple, il semble que votre code soit à peu près le même que:
for val in itertools.product(np.arange(0, 1, 0.1), repeat=6):
print F(val)
tous ces itérateurs produisent des tuples, pas des listes ou des tableaux numpy, donc si votre F est difficile à obtenir spécifiquement un tableau numpy, vous devrez accepter la charge supplémentaire de construire ou de nettoyer et remplir un à chaque étape.
numpy mise en œuvre doit être env. 2x la vitesse de la réponse donnée:
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
il semble que vous vouliez une grille pour évaluer votre fonction, dans ce cas vous pouvez utiliser numpy.ogrid (ouvert) ou numpy.mgrid (étoffé):
import numpy
my_grid = numpy.mgrid[[slice(0,1,0.1)]*6]
Vous pouvez faire quelque chose comme cela
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])
qui donne
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]])
Voici encore une autre façon, en utilisant pur NumPy, pas de récursion, pas de compréhension de liste, et pas d'explicite pour les boucles. Il est environ 20% plus lent que la réponse originale, et il est basé sur 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
par exemple,
x = np.arange(3)
a = cartesian(x, x, x, x, x)
print(a)
donne
[[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]]
pour une pure implémentation numpy du produit cartésien des tableaux 1D (ou des listes Python plates), il suffit d'utiliser meshgrid() , de rouler les axes avec transpose() , et de remodeler vers l'ouput désiré:
def cartprod(*arrays):
N = len(arrays)
return transpose(meshgrid(*arrays, indexing='ij'),
roll(arange(N + 1), -1)).reshape(-1, N)
Note Ceci a la convention du dernier axe changeant le plus rapidement ("style C" ou "ligne-majeur").
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]])
si vous voulez changer le premier axe le plus rapide ("FORTRAN style" ou "colonne-major"), il suffit de changer le order paramètre de reshape() comme ceci: reshape((-1, N), order='F')