Re: [R] sign(<permutation>) in R ?

From: Martin Maechler <maechler_at_stat.math.ethz.ch>
Date: Wed, 16 Apr 2008 18:51:12 +0200

As alluded, I've got several off-line replies and suggestions, notably Tobias Verbeke mentioning Matlab code {as an example of "pseudo code" :-) } by Burkhardt on the web,
and Barry Rowlinson providing code similar to Peter Dalgaard's, see below.
I've found myself python code on the web, unfortunately bugous (only for *some* permutation vectors).

Of all I've got, Peter's code was clearly the best, when n <- length(<permutation>) was not small. For that reason, I've ``speeded'' Peter's code, and ended in code that won hands down.

The simulation-testing of the diverse algorithms has been fun, not the least finding that one "public" algorithm was faulty; If I have a bit more time, I should write a small R news article on this.
Here's the current best code
{and of course, a C version of that will be yet an order of  magnitude faster} :

## MM: Peter's function is so compact and elegant, ## now try to optimize for speed() : cycLengths2 <- function(p) {

    n <- length(p)
    x <- integer(n)
    ii <- seq_len(n)
    for (i in ii) {

	z <- ii[!x][1] # index of first unmarked x[] entry
	if (is.na(z)) break
	repeat { ## mark x[] <- i  for those in i-th cycle
	    x[z] <- i
	    z <- p[z]
	    if (x[z]) break
	}

    }
## Now, table(x) gives the cycle lengths,
## where split(seq_len(n), x) would give the cycles list
    tabulate(x, i - 1L) ## is quite a bit faster than the equivalent
## table(x)

}

## MM: Here's the transformation from the cycle lengths' vector to sign(<p>): signCycLength <- function(clen) 1L - (sum(clen %% 2 == 0) %% 2)*2L

for the *sign* of a permutation 'p'
the above has to be called as signCycLength(cycLengths2(p)) e.g.,

   p <- sample(50)
   signCycLength(cycLengths2(p))

Thank you, again for all the feedback, and suggestions, Martin Maechler, ETH Zurich

>>>>>> Peter Dalgaard wrote yesterday

>> Martin Maechler wrote:
>> > Ok,
>> > thanks a lot, everyone!
>> >
>> > Yes, I should rather have started thinking a bit more myself,
>> > before going the easy route to R-help....
>> >
>> > Anyway, the most obvious algorithm,
>> > just putting things into place by swapping elements,
>> > and counting how many times you have to swap, is easy and
>> > quite efficient.
>> >
>> > I'll post R code later, being busy for the next few hours.
>> > Martin
>> >
>> >   
>> It is also quite easy to generate the cycles explicitly by a marking 
>> algorithm:
>> 
>> p <- sample(1:100)
>> x <- integer(100)
>> for (i in 1:100) {
>>     z <- which(!x)[1]
>>     if (is.na(z)) break
>>     repeat {
>>         x[z] <- i
>>         z <- p[z]
>>         if (x[z]) break
>>     }
>> }
>> table(x)
>> 
>> (the which(!x)[1] bit could be optimized somewhat if this had been C. 
>> Notice that the loop is essentially O(N) because every iteration marks 
>> one element of x)
>> 
>> 
>> >   
>> >>>>>> "MM" == Martin Maechler <maechler_at_stat.math.ethz.ch>
>> >>>>>>     on Tue, 15 Apr 2008 18:13:43 +0200 writes:
>> >>>>>>             
>> >
>> >     MM> I am looking for an algorithm (written in R (preferably) or C,
>> >     MM> but even pseudo-code in a text book maybe fine)
>> >     MM> to determine the sign of a permutation.
>> >
>> >     MM> What is that?  Well, a permutation is either even or odd, the
>> >     MM> sign is +1 or -1, respectively, see, e.g.,
>> >     MM> http://en.wikipedia.org/wiki/Signature_of_a_permutation
>> >     MM> which also says
>> >     >>> In practice, in order to determine whether a given permutation
>> >     >>> is even or odd, one writes the permutation as a product of
>> >     >>> disjoint cycles. The permutation is odd if and only if this
>> >     >>> factorization contains an odd number of even-length cycles.
>> >
>> >     MM> but I would not know how to algorithmically
>> >     MM> "write the permutation as a product of disjoint cycles"
>> >
>> >     MM> If you start looking at R code, 
>> >     MM> let's assume the permutation {\pi(i)}_{i=1..n} is simply given
>> >     MM> as the (integer) vector (\pi(1), \pi(2), ..., \pi(n))
>> >     MM> {or equivalently, a random permutation is simply found by  'sample(n)'}
>> >
>> >     MM> Thank you in advance for further pointers,
>> >     MM> or even working R code.
>> >
>> >     MM> Best regards,
>> >     MM> Martin Maechler, ETH Zurich
>> >

>> 
>> -- 
>>    O__  ---- Peter Dalgaard             Řster Farimagsgade 5, Entr.B
>>   c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
>>  (*) \(*) -- University of Copenhagen   Denmark      Ph:  (+45) 35327918
>> ~~~~~~~~~~ - (p.dalgaard_at_biostat.ku.dk)              FAX: (+45) 35327907

______________________________________________
R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Wed 16 Apr 2008 - 17:00:27 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Wed 16 Apr 2008 - 18:30:28 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive