Python - Случайный разрыв связей при выборе режима

Question

scipy.stats.mode прекрасно работает , но мне нужно разорвать модальные связи случайным образом.

import numpy as np
import scipy.stats as stats

a = np.array([[3, 3, 4], 
              [3, 1, 0], 
              [4, 5, 0], 
              [4, 3, 0]])

stats.mode(a, axis=0)

Out[37]: ModeResult(mode=array([[3, 3, 0]]), count=array([[2, 2, 3]]))

Для первого результата (столбец), scipy.stats.mode выбирает 3 среди связанных кандидаты 3 и 4, следующим образом :

Если существует более одного такого значения, возвращается только самое маленькое.

Таким образом, среди 3 и 4, он выбирает 3, потому что это самый маленький. Я хотел бы выбрать случайным образом 3 и 4, но scipy.stats.mode не возвращает достаточно информации, чтобы позволить мне это сделать. Есть ли хороший способ сделать это, используя numpy или достойную альтернативу?

Answer 1

Для более эффективного подхода вот альтернатива numba:

from numba import njit, int32

@njit
def mode_rand_ties(a):
    out = np.zeros(a.shape[1], dtype=int32)
    for col in range(a.shape[1]):
        z = np.zeros(a[:,col].max()+1, dtype=int32)
        for v in a[:,col]:
            z[v]+=1
        maxs = np.where(z == z.max())[0]
        out[col] = np.random.choice(maxs)
    return out

---

Там, где тестирование для массива выше, запустив несколько раз, мы увидим, что мы можем получить либо 3, либо 4 как режим первого столбца:

mode_rand_ties(a)
# array([4, 3, 0], dtype=int32)

mode_rand_ties(a)
# array([3, 3, 0], dtype=int32)

---

И, проверив производительность массива в форме (4000, 3), мы получаем, что это занимает всего около 40us:

x = np.concatenate([a]*1000, axis=0)
%timeit mode_rand_ties(x)
# 41.1 µs ± 13.2 µs per loop (mean ± std. dev. of 7 runs, 1 loop each)

Принимая во внимание, что с текущим решением:

%timeit mode_rand(x, axis=0)
# 388 µs ± 23.2 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

Answer 2

Я все еще отвечаю, если кто-нибудь придумает лучший способ, но это мое давнее временное решение, которое просто искажает исходный код scipy.stats.mode. Единственная существенная модификация - в for ind in inds l oop, где я использую np.where, чтобы вернуть все индексы с одинаковым числом максимальных значений, и я случайным образом выбираю индекс из этого.

from collections import namedtuple
ModeResult = namedtuple('ModeResult', ('mode', 'count'))
def mode_rand(a, axis):
    in_dims = list(range(a.ndim))
    a_view = np.transpose(a, in_dims[:axis] + in_dims[axis+1:] + [axis])

    inds = np.ndindex(a_view.shape[:-1])
    modes = np.empty(a_view.shape[:-1], dtype=a.dtype)
    counts = np.zeros(a_view.shape[:-1], dtype=np.int)

    for ind in inds:
        vals, cnts = np.unique(a_view[ind], return_counts=True)
        maxes = np.where(cnts == cnts.max())  # Here's the change
        modes[ind], counts[ind] = vals[np.random.choice(maxes[0])], cnts.max()

    newshape = list(a.shape)
    newshape[axis] = 1
    return ModeResult(modes.reshape(newshape), counts.reshape(newshape))

mode_rand(a, axis=0)