From r.hankin at noc.soton.ac.uk  Fri Jun  9 10:25:50 2006
From: r.hankin at noc.soton.ac.uk (Robin Hankin)
Date: Fri, 9 Jun 2006 09:25:50 +0100
Subject: [Rd] Idempotent apply
In-Reply-To: <f8e6ff050606090050t11b677aes463a877f7f61a272@mail.gmail.com>
References: <f8e6ff050606090050t11b677aes463a877f7f61a272@mail.gmail.com>
Message-ID: <54E36AE0-DE4E-4B96-8863-452DA75900EE@soc.soton.ac.uk>

Hi Hadley, Rteam,

I encountered similar issues trying to reproduce Matlab's "sort"
functionality.

[here  "a" is array(1:27, c(2,3,4))]

In matlab, one can type sort(a,i) to sort along dimension i.

One can reproduce this functionality thus:

asort <- function(a,i,FUN="sort"){
   j <- 1:length(dim(a))
   aperm(apply(a,j[-i],FUN),append(j[-1],1,i-1))
}


then  things like

asort(a,1)
asort(a,2)
asort(a,1,sample)

work as expected (at least by this recovering Matlab user)



best wishes



rksh





On 9 Jun 2006, at 08:50, hadley wickham wrote:

> Dear all,
>
> I have been working on an idempotent version of apply, such that
> applying a function f(x) = x (ie. force) returns the same array (or a
> permutation of it depending on the order of the margins):
>
> a <- array(1:27, c(2,3,4))
>
> all.equal(a, iapply(a, 1, force))
> all.equal(a, iapply(a, 1:2, force))
> all.equal(a, iapply(a, 1:3, force))
> all.equal(aperm(a, c(2,1,3)), iapply(a, 2, force))
> all.equal(aperm(a, c(3,1,2)), iapply(a, 3, force))
>
> I have also generalised apply so that the function can return an array
> of any dimension and it will be slotted in orthogonally to the
> existing dimensions:
>
> iapply(a, 1, min)
> iapply(a, 2, min)
> iapply(a, 3, min)
> iapply(a, 1, range)
> iapply(a, 2, range)
> iapply(a, 3, range)
>
> I have included my code below.  I'd be interested to get your  
> feedback on:
>
>  * whether you can find an edge case where it doesn't work
>
>  * how I can make the subsetting code more elegant - the current
> kludgework of do.call seems to be suboptimal, but I couldn't figure
> out a better way
>
>  * I was also suprised that a & b display differently in this example:
>
> a <- b <- as.array(1:3)
> dim(b) <- c(3,1)
>
> Any other comments are very much appreciated!
>
> Hadley
>
>
> iapply <- function(x, margins=1, fun, ..., REDUCE=TRUE) {
> 	if (!is.array(x)) x <- as.array(x)
> 	
> 	reorder <- c(margins, (1:length(dim(x)))[-margins])
>
> 	x <- aperm(x, (reorder))
> 	x <- compactify(x, length(margins))
>
> 	results <- lapply(x, fun, ...)
> 	dim(results) <- dim(x)
> 	
> 	results <- decompactify(results)
> 	if (REDUCE) reduce(results) else results
> }
> vdim <- function(x) if (is.vector(x)) length(x) else dim(x)
>
>
> # Compacts an array into a smaller array of lists containing the
> remaining dimensions
> compactify <- function(x, dims = 1) {
>
> 	d <- dim(x)
> 	ds <- seq(along=d)
> 	margins <- 1:dims
> 	
> 	subsets <- do.call(expand.grid, structure(lapply(d[margins], seq,
> from=1), names=margins))
> 	subsets[, ds[-margins]] <- TRUE
> 	
> 	res <- lapply(1:nrow(subsets), function(i) do.call("[",c(list(x),
> subsets[i, ], drop=TRUE)))
> 	dim(res) <- dim(x)[margins]
> 	
> 	res
> }
>
> # Inverse of compactity
> decompactify <- function(x) {
>
> 	subsets <- do.call(expand.grid, structure(lapply(dim(x), seq,  
> from=1)))
> 	subsets <- cbind(subsets, x=matrix(TRUE, ncol=length(vdim(x[[1]])),
> nrow=nrow(subsets)))
>
> 	y <- array(NA, c(vdim(x), vdim(x[[1]])))
> 	for(i in 1:length(x)) {
> 		y <- do.call("[<-",c(list(y), unname(subsets[i, ]), value = list(x 
> [[i]])))
> 	}
> 	y
> }	
> reduce <- function(x) {
> 		do.call("[", c(list(x), lapply(dim(x), function(x) if (x==1) 1 else
> T), drop=TRUE))
> }
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel

--
Robin Hankin
Uncertainty Analyst
National Oceanography Centre, Southampton
European Way, Southampton SO14 3ZH, UK
  tel  023-8059-7743



[Rd] Idempotent apply
Robin Hankin