[Numpy-discussion] Permutations in Simulations`

[Mark Miller]

Out of curiosity, why wouldn't numpy.apply_along_axis be a reasonable
approach here.  Even more curious:  why is it slower than the original
explicit loop?

Learning,

-Mark

import numpy as np
import timeit

def baseshuffle(nx, ny):
    x = np.arange(nx)
    res = np.zeros((nx,ny),int)
    for sim in range(ny):
        res[:,sim] = np.random.permutation(x)
    return res

def axisshuffle(nx,ny):
    x=np.arange(nx).reshape(nx,1)
    res=np.zeros((nx,ny),int)
    res[:] = x
    x1 = np.apply_along_axis(np.random.permutation,0,res)
    return x1

nx = 10
ny = 100

t=timeit.Timer("baseshuffle(nx, ny)", "from __main__ import *")
print t.timeit(100)
t=timeit.Timer("axisshuffle(nx,ny)", "from __main__ import *")
print t.timeit(100)


>>>
0.111320049056
0.297792159759
>>>












On Tue, Feb 10, 2009 at 12:58 PM, Keith Goodman <kwgoodman at gmail.com> wrote:
> On Tue, Feb 10, 2009 at 12:41 PM, Keith Goodman <kwgoodman at gmail.com> wrote:
>> On Tue, Feb 10, 2009 at 12:28 PM, Keith Goodman <kwgoodman at gmail.com> wrote:
>>> On Tue, Feb 10, 2009 at 12:18 PM, Keith Goodman <kwgoodman at gmail.com> wrote:
>>>> On Tue, Feb 10, 2009 at 11:29 AM, Mark Janikas <mjanikas at esri.com> wrote:
>>>>> I want to create an array that contains a column of permutations for each
>>>>> simulation:
>>>>>
>>>>> import numpy as NUM
>>>>>
>>>>> import numpy.random as RAND
>>>>>
>>>>> x = NUM.arange(4.)
>>>>>
>>>>> res = NUM.zeros((4,100))
>>>>>
>>>>>
>>>>> for sim in range(100):
>>>>>
>>>>> res[:,sim] = RAND.permutation(x)
>>>>>
>>>>>
>>>>> Is there a way to do this without a loop?  Thanks so much ahead of time…
>>>>
>>>> Does this work? Might not be faster but it does avoid the loop.
>>>>
>>>> import numpy as np
>>>>
>>>> def weirdshuffle(nx, ny):
>>>>    x = np.ones((nx,ny)).cumsum(0, dtype=np.int) - 1
>>>>    yidx = np.ones((nx,ny)).cumsum(1, dtype=np.int) - 1
>>>>    xidx = np.random.rand(nx,ny).argsort(0).argsort(0)
>>>>    return x[xidx, yidx]
>>>
>>> Hey, it is faster for nx=4, ny=100
>>>
>>> def baseshuffle(nx, ny):
>>>    x = np.arange(nx)
>>>    res = np.zeros((nx,ny))
>>>    for sim in range(ny):
>>>        res[:,sim] = np.random.permutation(x)
>>>    return res
>>>
>>>>> timeit baseshuffle(4,100)
>>> 1000 loops, best of 3: 1.11 ms per loop
>>>>> timeit weirdshuffle(4,100)
>>> 10000 loops, best of 3: 127 µs per loop
>>>
>>> OK, who can cut that time in half? My first try looks clunky.
>>
>> This is a little faster:
>>
>> def weirdshuffle2(nx, ny):
>>    one = np.ones((nx,ny), dtype=np.int)
>>    x = one.cumsum(0)
>>    x -= 1
>>    yidx = one.cumsum(1)
>>    yidx -= 1
>>    xidx = np.random.random_sample((nx,ny)).argsort(0).argsort(0)
>>    return x[xidx, yidx]
>>
>>>> timeit weirdshuffle(4,100)
>> 10000 loops, best of 3: 129 µs per loop
>>>> timeit weirdshuffle2(4,100)
>> 10000 loops, best of 3: 106 µs per loop
>
> Sorry for all the mail.
>
> def weirdshuffle3(nx, ny):
>    return np.random.random_sample((nx,ny)).argsort(0).argsort(0)
>
>>> timeit weirdshuffle(4,100)
> 10000 loops, best of 3: 128 µs per loop
>>> timeit weirdshuffle3(4,100)
> 10000 loops, best of 3: 37.5 µs per loop
> _______________________________________________
> Numpy-discussion mailing list
> Numpy-discussion at scipy.org
> http://projects.scipy.org/mailman/listinfo/numpy-discussion
>