Services & ResourcesWolfram ForumsMathGroup Archive

[Date Index] [Thread Index] [Author Index]

RE: generate random permutation

• To: mathgroup at smc.vnet.net
• Subject: [mg40319] RE: [mg40287] generate random permutation
• From: "Wolf, Hartmut" <Hartmut.Wolf at t-systems.com>
• Date: Tue, 1 Apr 2003 04:51:05 -0500 (EST)
• Sender: owner-wri-mathgroup at wolfram.com

>-----Original Message-----
>From: Umlaut [mailto:r.dorne at ntlworld.com]
To: mathgroup at smc.vnet.net
>Sent: Sunday, March 30, 2003 11:08 AM
>To: mathgroup at smc.vnet.net
>Subject: [mg40319] [mg40287] generate random permutation
>
>
>
>I am a researcher in Artificial Intelligence and I have a 
>question about
>generating random permutation of a set of numbers.
> 
>More formal definition:
>Given a set of N numbers S = {n1, ., nn}
> 
>how to efficiently and accurately generate P(S) = random_permutation(S)
>Ex: S = {3,4,12}, P(S) should have a probability equals to 
>1/3! = 1/6 to
>return one of the following permutations:
>{(3,4,12)(3,12,4)(4,3,12)(4,12,3)(12,3,4)(12,4,3)}
>and a complexity equals to O(N).
> 
>and
> 
>how to efficiently generate an iterator on the numbers which will be
>returned by P(S) ?
>Ex: S = {3,4,12} P(S) can return any of the permutations described
>below, but what we want here is to get an iterator on the 
>future numbers 
>without creating the full permutation 
> 
> 
>Do you know any algorithm performing this task which is very efficient
>in terms of quality (each permutation has the same probability to be 
>generated) and in terms of complexity (minimal number of operations
>is necessary to perform this operation n2, n, log n, etc.).
> 
>If you have any clue or web links where I can find this information, 
>this would be very much appreciated.
> 
> 
>

The problem with that question are the conflicting goals of correctness,
performance and asymptotc behavior. The right answer depends on what you
want to do.


In Mathematica there are two packages:

In[1]:= << DiscreteMath`Permutations`

In[2]:= RandomPermutation[13070]; // Timing
Out[2]= {0.501 Second, Null}

In[3]:= Quit[]



And also 

In[1]:= << DiscreteMath`Combinatorica`

In[2]:= RandomPermutation1[13070]; // Timing
Out[2]= {0.42 Second, Null}

In[3]:= RandomPermutation2[13070]; // Timing
Out[3]= {1.132 Second, Null}

In[6]:= Quit[]



The fastest algorithm, I know of is this:

In[1]:=
permutation[n_Integer?NonNegative] := Ordering[Array[Random[] &, {n}]]

In[2]:= permutation[13070]; // Timing
Out[2]= {0.04 Second, Null}
In[3]:= Quit[]



Depending on your problem, you may be happy with this. It is however
somewhat flawed, as well as RandomPermutation and RandomPermutation1: If the
same random reals are hit twice (improbable, though possible) the Sort or
Ordering will not give these two in random order. Now it depends on what is
more important: this rather seldom incidence or performance. You may try to
eliminate the problem, by checking e.g. with Split. A penalty has to be
paid, of course:

In[1]:= 
permutation2[n_Integer?NonNegative] := 
  Block[{ra = Array[Random[] &, {n}], ord},
    If[Times @@ Length /@ Split[ra[[ord = Ordering[ra]]]] === 1, ord, 
      permutation2[n]]]

In[2]:= permutation2[13070]; // Timing
Out[2]= {0.18 Second, Null}
In[3]:= Quit[]

(Of course this discussion also questions the quality of the random number
generator behind, such that this "improvement" may become nonsense for other
reasons.)

The asymptotic behaviour of these algorithms is O[n log n], i.e. good for
all practical purposes.



Another point is whether is is senseful to have an iterator to spit out the
elements of the random permutation one by one.  Im Mathematica I'd say no,
list operation is faster, however (and depending on your real problem) in
C++ things might become different. Here is an iterator for Mathematica:

 
In[1]:=
makeRPI /: Set[lhs_, makeRPI[n_]] := (lhs = makeRPI0[n];)

In[2]:= makeRPI0[n_Integer?Positive] :=
  Module[{r = Range[n], nr = n, x},
    If[nr === 0, {},
        {x, r} = Through[{Part, Delete}[r, Random[Integer, {1, nr--}]]]; x]
&]

In[3]:= it = makeRPI[5]

In[4]:= Table[it[], {6}]
Out[4]= {1, 5, 4, 2, 3, {}}

In[5]:=
(it2 = makeRPI[13070]; Table[it2[], {13070}];) // Timing
Out[5]= {68.819 Second, Null}
In[6]:= Quit[]

This algorithm for each step draws at random from a bag of Integers. The
Problem with this one is the Delete-ing of the integer drawn: this makes its
behaviour as O[n^2].



Instead of deleting, we just loose the integer drawn out of sight from the
window of interest, our bag now just being the rest of an initial list of
integers to be mapped over (the integer not drawn is moved into the position
of the integer drawn):

In[1]:=
randperm2[n_] := 
  Module[{deck = Range[n], s, j}, 
    MapIndexed[(s = deck[[j = #1 + First[#2] - 1]]; 
          deck[[j]] = deck[[First[#2]]]; s) &, 
      Reverse[Table[Random[Integer, {1, i}], {i, n}]]]]

In[2]:= randperm2[5]
Out[2]= {2, 5, 3, 4, 1}

In[3]:= randperm2[13070]; // Timing
Out[3]= {1.252 Second, Null}
In[4]:= Quit[]



>From this it's not a far way to the "shuffle algorithm" presented to this
community by Daniel Lichtblau (search the archive), which can be compiled: 

In[1]:=
randperm = 
    Compile[{{n, _Integer}}, 
      Module[ {deck = Range[n], newj},  
        Do[ newj = Random[Integer, {j, n}]; 
          deck[[{j, newj}]] = deck[[{newj, j}]], {j, n - 1}]; deck ]] ;

In[2]:= randperm[5]
Out[2]= {4, 3, 2, 1, 5}

In[3]:= randperm[13070]; // Timing
Out[3]= {0.21 Second, Null}


(RandomPermutation2 from the Package DiscreteMath`Combinatorica` is quite
similar, however not coded as effectively and not compiled.)


--
Hartmut Wolf