[R] plyr: a*ply with functions that return matrices-- possible bug in aaply?
    Michael Friendly 
    friendly at yorku.ca
       
    Sun Oct  3 22:04:44 CEST 2010
    
    
  
  I have an application where I have a function to calculate results for 
a 2-way table or matrix, which
returns a matrix with one less row and column. To keep this short, the 
function below captures the structure:
fun2way <- function(f){
     if (!length(dim(f)) ==2) stop("only for 2-way arrays")
     R <- dim(f)[1]
     C <- dim(f)[2]
     f[1:(R-1), 1:(C-1)]
}
Now, I want to extend this to higher-way arrays, using apply-like 
methods over the strata (all but the first two dimensions),
and returning an array in which the last dimensions correspond to 
strata.  That is, I want to define something like the
following using an a*ply method, but aaply gives a result in which the 
applied .margin(s) do not appear last in the
result, contrary to the documentation for ?aaply.  I think this is a 
bug, either in the function or the documentation,
but perhaps there's something I misunderstand for this case.
fun <- function(f, stratum=NULL) {
     L <- length(dim(f))
   if (L > 2 && is.null(stratum))
         stratum <- 3:L
   if (is.null(stratum)) {
     result <- fun2way(f)
   }
   else {
         require(plyr)
     result <- aaply(f, stratum, fun2way)  ## order of dimensions 
screwed up!
     }
   result
}
For example, by hand (or with a loop) I can calculate the
pieces and combine them as I want using abind():
 > # apply separately to strata
 > t1<-fun2way(HairEyeColor[,,1])
 > t2<-fun2way(HairEyeColor[,,2])
 >
 > library(abind)
 > abind(t1, t2, along=3)
, , 1
       Brown Blue Hazel
Black    32   11    10
Brown    53   50    25
Red      10   10     7
, , 2
       Brown Blue Hazel
Black    36    9     5
Brown    66   34    29
Red      16    7     7
alply() gives me what I want, but with the strata as list elements, 
rather than an array
 > library(plyr)
 > # strata define separate list elements
 > alply(HairEyeColor, 3, fun2way)
$`1`
        Eye
Hair    Brown Blue Hazel
   Black    32   11    10
   Brown    53   50    25
   Red      10   10     7
$`2`
        Eye
Hair    Brown Blue Hazel
   Black    36    9     5
   Brown    66   34    29
   Red      16    7     7
attr(,"split_type")
[1] "array"
attr(,"split_labels")
      Sex
1   Male
2 Female
 >
However, with aaply(), dim[3] ends up as first dimension, not last
 > # dim[3] ends up as first dimension, not last
 > aaply(HairEyeColor, 3, fun2way)
, , Eye = Brown
         Hair
Sex      Black Brown Red
   Female    36    66  16
   Male      32    53  10
, , Eye = Blue
         Hair
Sex      Black Brown Red
   Female     9    34   7
   Male      11    50  10
, , Eye = Hazel
         Hair
Sex      Black Brown Red
   Female     5    29   7
   Male      10    25   7
 > str(aaply(as.array(HairEyeColor), 3, fun2way))
  num [1:2, 1:3, 1:3] 36 32 66 53 16 10 9 11 34 50 ...
  - attr(*, "dimnames")=List of 3
   ..$ Sex : chr [1:2] "Female" "Male"
   ..$ Hair: chr [1:3] "Black" "Brown" "Red"
   ..$ Eye : chr [1:3] "Brown" "Blue" "Hazel"
 >
 > ## aaply should return this....
 > aperm(aaply(HairEyeColor, 3, fun2way), c(2,3,1))
, , Sex = Female
        Eye
Hair    Brown Blue Hazel
   Black    36    9     5
   Brown    66   34    29
   Red      16    7     7
, , Sex = Male
        Eye
Hair    Brown Blue Hazel
   Black    32   11    10
   Brown    53   50    25
   Red      10   10     7
 >
On the other hand, aaply() does work as I expect, with an array of size 
2 x C x strata
 > library(vcd)
 > fun2way(Employment[,,1])
<1Mo  1-3Mo 3-12Mo  1-2Yr  2-5Yr
      8     35     70     62     56
 > fun2way(Employment[,,2])
<1Mo  1-3Mo 3-12Mo  1-2Yr  2-5Yr
     40     85    181     85    118
 >
 > aaply(Employment, 3, fun2way)
LayoffCause <1Mo 1-3Mo 3-12Mo 1-2Yr 2-5Yr
    Closure     8    35     70    62    56
    Replaced   40    85    181    85   118
-- 
Michael Friendly     Email: friendly AT yorku DOT ca
Professor, Psychology Dept.
York University      Voice: 416 736-5115 x66249 Fax: 416 736-5814
4700 Keele Street    Web:   http://www.datavis.ca
Toronto, ONT  M3J 1P3 CANADA
    
    
More information about the R-help
mailing list