Python - Randomly breaking ties when choosing a mode

Question

scipy.stats.mode works great, but I need to break modal ties at random.

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]]))

For the first result (column), scipy.stats.mode chooses 3 among the tied candidates 3 and 4, as follows:

If there is more than one such value, only the smallest is returned.

So among 3 and 4, it picks 3 because it's the smallest. I'd like to randomly choose among 3 and 4, but scipy.stats.mode doesn't bring back enough information to allow me to do that. Is there a good way to do this using numpy or a decent alternative?

Answer 1

For a performant approach, here's a numba alternative:

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

---

Where testing for the array above, by running more than once we see that we can get either 3 or 4 as the mode of the first column:

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

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

---

And by checking the performance on a (4000, 3) shaped array we get that it takes only about 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)

Whereas with the current solution:

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