Szeregowanie tablic wielowymiarowych na miejscu

Question

Próbuję zaimplementować procedurę losowania bezpieczną dla NaN w języku Cython, która może tasować wzdłuż kilku osi wielowymiarowej macierzy o dowolnym wymiarze.

W prostym przypadku matrycy 1D, może po prostu losowo na wszystkie współczynniki o braku NaN wartości za pomocą algorytmu Fisher Yates:

def shuffle1D(np.ndarray[double, ndim=1] x): 
    cdef np.ndarray[long, ndim=1] idx = np.where(~np.isnan(x))[0] 
    cdef unsigned int i,j,n,m 

    randint = np.random.randint 
    for i in xrange(len(idx)-1, 0, -1): 
     j = randint(i+1) 
     n,m = idx[i], idx[j] 
     x[n], x[m] = x[m], x[n] 

ja Aby rozszerzyć ten algorytm obsługiwać duże wielowymiarowy tablice bez zmiany kształtu (które uruchamiają kopię dla bardziej skomplikowanych przypadków, których tutaj nie uwzględniono). W tym celu musiałbym pozbyć się stałego wymiaru wejściowego, co nie wydaje się możliwe w przypadku numpy array ani memoryviews w Cython. Czy jest w pobliżu praca?

Wielkie dzięki z góry!

Answer 1

Dzięki komentarzach @Veedrac ta odpowiedź wykorzystuje więcej możliwości Cython.

• tablicy wskaźnik przechowuje adres pamięci wartościami wzdłuż axis
• Twój algorytm jest używany z modyfikacją that checks for nan values, zapobiegając ich być klasyfikowane
• To nie utworzy kopię C zamówionych tablic. W przypadku uporządkowanych tablic Fortran, polecenie ravel() zwróci kopię. To można poprawić poprzez tworzenie kolejną tablicę podwójnych wskaźników do przewozu wartości x, prawdopodobnie z jakiejś kary cache ...

Ten kod jest co najmniej jeden rząd wielkości szybciej niż inne oparte na plasterki.

from libc.stdlib cimport malloc, free 

cimport numpy as np 
import numpy as np 
from numpy.random import randint 

cdef extern from "numpy/npy_math.h": 
    bint npy_isnan(double x) 

def shuffleND(x, int axis=-1): 
    cdef np.ndarray[double, ndim=1] v # view of x 
    cdef np.ndarray[int, ndim=1] strides 
    cdef int i, j 
    cdef int num_axis, pos, stride 
    cdef double tmp 
    cdef double **v_axis 

    if axis==-1: 
     axis = x.ndim-1 

    shape = list(x.shape) 
    num_axis = shape.pop(axis) 

    v_axis = <double **>malloc(num_axis*sizeof(double *)) 
    for i in range(num_axis): 
     v_axis[i] = <double *>malloc(1*sizeof(double)) 

    try: 
     tmp_strides = [s//x.itemsize for s in x.strides] 
     stride = tmp_strides.pop(axis) 
     strides = np.array(tmp_strides, dtype=np.int32) 
     v = x.ravel() 
     for indices in np.ndindex(*shape): 
      pos = (strides*indices).sum() 
      for i in range(num_axis): 
       v_axis[i] = &v[pos + i*stride] 
      for i in range(num_axis-1, 0, -1): 
       j = randint(i+1) 
       if npy_isnan(v_axis[i][0]) or npy_isnan(v_axis[j][0]): 
        continue 
       tmp = v_axis[i][0] 
       v_axis[i][0] = v_axis[j][0] 
       v_axis[j][0] = tmp 
    finally: 
     free(v_axis) 

    return x 

Answer 2

Poniższy algorytm jest oparty na plasterkach, gdzie nie jest wykonywana żadna kopia i powinna działać dla każdego np.ndarray. Główne etapy to:

• np.ndindex() służy do uruchamiania throught różnych indeksów wielowymiarowych, z wyjątkiem jednego należącej do osi chcesz przetasować
• shuffle już opracowaną przez Ciebie dla przypadku 1-D jest stosowana .

Kod:

def shuffleND(np.ndarray x, axis=-1): 
    cdef np.ndarray[long long, ndim=1] idx 
    cdef unsigned int i, j, n, m 
    if axis==-1: 
     axis = x.ndim-1 
    all_shape = list(np.shape(x)) 
    shape = all_shape[:] 
    shape.pop(axis) 
    for slices in np.ndindex(*shape): 
     slices = list(slices) 
     axis_slice = slices[:] 
     axis_slice.insert(axis, slice(None)) 
     idx = np.where(~np.isnan(x[tuple(axis_slice)]))[0] 
     for i in range(idx.shape[0]-1, 0, -1): 
      j = randint(i+1) 
      n, m = idx[i], idx[j] 
      slice1 = slices[:] 
      slice1.insert(axis, n) 
      slice2 = slices[:] 
      slice2.insert(axis, m) 
      slice1 = tuple(slice1) 
      slice2 = tuple(slice2) 
      x[slice1], x[slice2] = x[slice2], x[slice1] 
    return x