Stochastically rounding a float to an integer

Question

I want a function (using Python 3.6+ if it's relevant) that will stochastically round a floating-point number to an integer in the following manner:

Given a real number, x, let a = floor(x) and let b = ceil(x). Then, write a function s_int() that will return a with a probability of b - x and return b with a probability of x - a.

For example, s_int(14.8) should return 14 20% of the time, and return 15 for the remaining 80% of the time.

---

Here is my attempt:

import math
from random import random

def s_int(x):
    a = math.floor(x)
    return a + ((x - a) > random())

It appears to work for all cases I can think of:

In [2]: Counter(s_int(14.7) for _ in range(1000000))
Out[2]: Counter({14: 300510, 15: 699490})

In [3]: Counter(s_int(-14.7) for _ in range(1000000))
Out[3]: Counter({-15: 700133, -14: 299867})

In [4]: Counter(s_int(14) for _ in range(1000000))
Out[4]: Counter({14: 1000000})

In [5]: Counter(s_int(-14) for _ in range(1000000))
Out[5]: Counter({-14: 1000000})

In [6]: Counter(s_int(0) for _ in range(1000000))
Out[6]: Counter({0: 1000000})

---

Here are my questions:

1. Are there any edge cases I didn't consider that this function doesn't work for?

2. Are there other simpler or more elegant solutions?

3. Can this be made to run faster?

Answer 1

1. I think there are no edge cases
2. Try using numpy

import numpy as np
def s_int(x):
    a = np.floor(x)
    b = a + 1
    return (np.random.choice([a, b], p=[b - x, x - a]))

3. I think that it's an O(1) operation. Does not go faster.