0

EDIT: I have simplified the example and the functions following the comment of @jdorbes. Hope this improves the post.

I have a series of N-dimensional array that I need to "rotate" along their first dimension (e.g. in a 2d matrix, switching between each other the first and the second half of the columns). For example considering a 3D case I can define my array as:

field <- array(1:32, c(4, 4, 2))

I developed a very basic and dimension-dependent function that achieves the goal. The same operation can be done of course also with head and tail but I found out that this is the faster way:

rotate.test <- function(field) {
  xxx <- field
  dims <- length(dim(xxx))
  if (dims == 2) { # for x,y data
      ll <- length(field[, 1])
      field[(ll * 0.5 + 1  ):ll, ] <- xxx[1:(ll * 0.5), ]
      field[1:(ll * 0.5), ] <- xxx[(ll * 0.5  + 1):ll, ]
  }
  if (dims == 3) { # for x,y,t data
      ll <- length(field[, 1, 1]) 
      field[(ll * 0.5 + 1 ):ll, , ] <- xxx[1:(ll * 0.5), , ]
      field[1:(ll * 0.5), , ] <- xxx[(ll * 0.5 + 1):ll, , ]
  } 
    return(field)
}

The result reads as:

> rotate.test(field)
, , 1

     [,1] [,2] [,3] [,4]
[1,]    3    7   11   15
[2,]    4    8   12   16
[3,]    1    5    9   13
[4,]    2    6   10   14

, , 2

     [,1] [,2] [,3] [,4]
[1,]   19   23   27   31
[2,]   20   24   28   32
[3,]   17   21   25   29
[4,]   18   22   26   30

This can be further generalized for N-dimensional data, but I guess you see that it is very clumsy. Following this post Select along one of n dimensions in array I found out that I can generalize it introducing a very useful call to the do.call function:

array_indexing <- function(field, dim, value, drop = FALSE) {
  indices <- rep(list(bquote()), length(dim(field)))
  indices[[dim]] <- value
  out <- do.call("[",c(list(field), indices, list(drop = drop)))
  return(out)
}

Then using abind

rotate.evo <- function(field) {
 require("abind")
 xdim <- 1
 ll <- dim(field)[xdim]
 firstchunk <- array_indexing(field, xdim, (ll * 0.5 + 1):ll)
 secondchunk <- array_indexing(field, xdim, 1:(ll * 0.5))
 out <- abind(firstchunk, secondchunk, along = 1)
 return(out)
}

However, the benchmarking (done with rbenchmark) is awful likely due to the call to abind instead of the simpler index replacement.

                test replications elapsed relative user.self sys.self
2  rotate.evo(field)         1000   0.547     6.36     0.542    0.005
1 rotate.test(field)         1000   0.086     1.00     0.069    0.016
  user.child sys.child
2          0         0
1          0         0

Here it comes my question: is there any way to generalize the index replacement as done by the array_indexing function? In other words, is there any way to perform this switching in a fast and robust way without writing ad-hoc cases for each dimension? What I am trying to find out is exists in R a method to generalize the replacement of some array index along a defined dimension.

I know this may sound as a weird request, but it will be great if I can avoid to load extra packages. This code is included in a online tool that should be lighter as possible (actually, I am trying to get rid of abind too!)

Thanks for any hint, and please ask if any part of my explanation is incomplete. Best, Paolo

  • It's very difficult to tell from your post what the expected effect of your custom function is meant to be. Could you provide an example output using a simpler input, perhaps an array such as `array(1:32, c(4, 4, 2))`? – jdobres Jan 02 '19 at 17:53
  • @jdobres I tried to simplify the whole story. I am trying to find out if exists in R a method to generalize the replacement of a part of an array through indexing along a defined dimension without writing a dimension-dependent code. – Paolo Davini Jan 03 '19 at 07:44

0 Answers0