sampling integers uniformly efficiently in python using numpy/scipy
Question: I have a problem where depending on the result of a random coin flip, I have
to sample a random starting position from a string. If the sampling of this
random position is uniform over the string, I thought of two approaches to do
it: one using multinomial from numpy.random, the other using the simple
randint function of Python standard lib. I tested this as follows:

    
    
    from numpy import *
    from numpy.random import multinomial
    from random import randint
    import time
    
    def use_multinomial(length, num_points):
        probs = ones(length)/float(length)
        for n in range(num_points):
        result = multinomial(1, probs)
    
    def use_rand(length, num_points):
        for n in range(num_points):
        rand(1, length)
    
    def main():
        length = 1700
        num_points = 50000
    
        t1 = time.time()
        use_multinomial(length, num_points)
        t2 = time.time()
        print "Multinomial took: %s seconds" %(t2 - t1)
    
        t1 = time.time()
        use_rand(length, num_points)
        t2 = time.time()
        print "Rand took: %s seconds" %(t2 - t1)    
    
    if __name__ == '__main__':
        main()
    

The output is:

Multinomial took: 6.58072400093 seconds Rand took: 2.35189199448 seconds

it seems like randint is faster, but it still seems very slow to me. Is there
a vectorized way to get this to be much faster, using numpy or scipy?

thanks.


Answer: I changed your code to actually return values (and used `randint` instead of
`rand` \- isn't that what you meant?) like this...

    
    
    def use_multinomial(length, num_points):
        probs = ones(length)/float(length)
        return multinomial(1, probs, num_points)
    
    def use_rand(length, num_points):
        return [randint(1,length) for _ in range(num_points)]
    

Then I tried my own version, using `numpy.random.randint` to generate a numpy
array of random points on the string:

    
    
    def use_np_randint(length, num_point):
        return nprandint(1, length, num_points)
    

The results:

    
    
    Multinomial took: 13.6279997826 seconds
    Rand took: 0.185000181198 seconds
    NP randint took: 0.00100016593933 seconds
    

Multinomial is obviously really slow comparitively, but is that even what you
want? I thought you said you wanted a uniform distribution? Using numpy's
randint is clearly the fastest of the bunch.

