## Paul E. Johnson
## 2013-07-10

## I like 2 dimensional arrays: matrices. I like data frames.  Two
## separate clients have asked for array help this week. If you saw
## data-array-1.R, you see I made some progress. Maybe then I'll be
## more ready to help when a problem like that comes along again if I
## work out some more array examples.

## But I still hate it. I'm scouring Internet for examples to chew on.

##
## Example 1. a 3x3x9 matrix that needs to be grouped by values of
## a factor.

## http://stackoverflow.com/questions/16135877/applying-a-function-to-a-multidimensional-array-with-grouping-variable

V <- 1:81
dim(V) <- c(3, 3, 9) ## fast way to create array V
group <- c('a','a','a','b','b','b','c','c','c')
## "Given that the grouping variable has 3 levels (a, b and c),
## the result (out) I'm looking for is an array of dimension 3x3x3."

## Here's a solution that depends on the package abind

library(abind)

out1 <- apply(V[ , , c(1:3)], c(1, 2), sum)
out2 <- apply(V[ , , c(4:6)], c(1, 2), sum)
out3 <- apply(V[ , , c(7:9)], c(1, 2), sum)

out <- abind(out1, out2, out3, along = 3)

## User asks for a general R solution, not depending on abind.

## I have no idea what this code means.
out <- apply(V, c(1, 2), by, group, sum)

## The group and sum arguments aren't named, I can't tell what they are
## used for. So I had to go run by a whole lot of times, to figure
## out what it is doing.

## Aha. syntax: by(x, INDICES = group, FUN = sum)
## and x is a sub array extracted from V. Headache.

## Take one example from the array, choosing this piece:

myx1 <- V[1, 1, ]
## Apply by to just that one bit
by(myx1, INDICES = group, FUN = sum)

## group: a
## [1] 30
## --------------------------------------------------------
## group: b
## [1] 111
## --------------------------------------------------------
## group: c
## [1] 192

## Take another piece
myx2 <- V[1, 2, ]
by(myx2, INDICES = group, FUN = sum)

## group: a
## [1] 39
## --------------------------------------------------------
## group: b
## [1] 120
## --------------------------------------------------------
## group: c
## [1] 201

## so the apply puzzle clear up. It cycles over

## This uses aperm to take same and re-organize the dimensions.
out <- aperm(apply(V, c(1, 2), by, group, sum), c(2, 3, 1))

## OK, so I think I understand that one.

## Example 2. This is about "compressing" 2 "sections" into 1.
## I think it helps to visualize it with named dimensions and columns and rows
V <- array(1:27, c(3, 3, 3), dimnames = list("myDim1" = c("A","B","C"),
"myDim2" = c("G","H","I"),
"myDim3" = c("X","Y","Z")))
V

## , , myDim3 = X
##
##       myDim2
## myDim1 G H I
##      A 1 4 7
##      B 2 5 8
##      C 3 6 9
##
## , , myDim3 = Y
##
##       myDim2
## myDim1  G  H  I
##      A 10 13 16
##      B 11 14 17
##      C 12 15 18
##
## , , myDim3 = Z
##
##       myDim2
## myDim1  G  H  I
##      A 19 22 25
##      B 20 23 26
##      C 21 24 27

## Names really help me see what is up here
## Now you can ask for pieces by name!
V["A","G", ]

##  X  Y  Z
##  1 10 19

V["A", , ]

##       myDim3
## myDim2 X  Y  Z
##      G 1 10 19
##      H 4 13 22
##      I 7 16 25

## Suppose you want to aggregate values "X" and "Y" by summing.
V[ , , "X"]

##       myDim2
## myDim1 G H I
##      A 1 4 7
##      B 2 5 8
##      C 3 6 9

V[ , , "Y"]

##       myDim2
## myDim1  G  H  I
##      A 10 13 16
##      B 11 14 17
##      C 12 15 18

V[ , , "X"] + V[ , , "Y"]

##       myDim2
## myDim1  G  H  I
##      A 11 17 23
##      B 13 19 25
##      C 15 21 27

## Now, replace the bit that was V[ , , c("X","Y")] with that new thing.
## Appears this works, but seems dumb to me.
V[ , , "X"] <- V[ , , "X"] + V[ , , "Y"]

Vnew <- V[ , , c("X","Z")]
Vnames <- dimnames(Vnew)
Vnames$myDim3[1] <- "XandY" dimnames(Vnew) <- Vnames Vnew  ## , , myDim3 = XandY ## ## myDim2 ## myDim1 G H I ## A 11 17 23 ## B 13 19 25 ## C 15 21 27 ## ## , , myDim3 = Z ## ## myDim2 ## myDim1 G H I ## A 19 22 25 ## B 20 23 26 ## C 21 24 27  ## There's got to be a better way to collapse 2 slices ## Example 4. Take 3, but with a missing. ## Now, a relevant problem. What if there is an NA in one spot? V["A", "G", "X"] <- NA V  ## , , myDim3 = X ## ## myDim2 ## myDim1 G H I ## A NA 17 23 ## B 13 19 25 ## C 15 21 27 ## ## , , myDim3 = Y ## ## myDim2 ## myDim1 G H I ## A 10 13 16 ## B 11 14 17 ## C 12 15 18 ## ## , , myDim3 = Z ## ## myDim2 ## myDim1 G H I ## A 19 22 25 ## B 20 23 26 ## C 21 24 27  V[ , , "X"]  ## myDim2 ## myDim1 G H I ## A NA 17 23 ## B 13 19 25 ## C 15 21 27  V[ , , "Y"]  ## myDim2 ## myDim1 G H I ## A 10 13 16 ## B 11 14 17 ## C 12 15 18  V[ , , "X"] + V[ , , "Y"]  ## myDim2 ## myDim1 G H I ## A NA 30 39 ## B 24 33 42 ## C 27 36 45  ## If you want a missing in there, OK. But what if you want the ## NA to be treated as zero. Harder. ## Put the "Z" part out of the way. Here's one way V[ , , -3]  ## , , myDim3 = X ## ## myDim2 ## myDim1 G H I ## A NA 17 23 ## B 13 19 25 ## C 15 21 27 ## ## , , myDim3 = Y ## ## myDim2 ## myDim1 G H I ## A 10 13 16 ## B 11 14 17 ## C 12 15 18  ## That's risky, 3 might be the wrong number. Can choose by name. V[ , , -which(dimnames(V)[[3]] == "Z")]  ## , , myDim3 = X ## ## myDim2 ## myDim1 G H I ## A NA 17 23 ## B 13 19 25 ## C 15 21 27 ## ## , , myDim3 = Y ## ## myDim2 ## myDim1 G H I ## A 10 13 16 ## B 11 14 17 ## C 12 15 18  VxAndy <- apply(V[ , , -which(dimnames(V)[[3]] == "Z")], MARGIN = c(1,2), FUN = sum, na.rm = T) VxAndy  ## myDim2 ## myDim1 G H I ## A 10 30 39 ## B 24 33 42 ## C 27 36 45  Vnew2 <- V[ , , c("X", "Z")] Vnew2[ , , "X"] <- VxAndy Vnames <- dimnames(Vnew2) Vnames$myDim3[1] <- "XandY"
dimnames(Vnew2) <- Vnames
Vnew2

## , , myDim3 = XandY
##
##       myDim2
## myDim1  G  H  I
##      A 10 30 39
##      B 24 33 42
##      C 27 36 45
##
## , , myDim3 = Z
##
##       myDim2
## myDim1  G  H  I
##      A 19 22 25
##      B 20 23 26
##      C 21 24 27

## Remember where we started.
Vnew

## , , myDim3 = XandY
##
##       myDim2
## myDim1  G  H  I
##      A 11 17 23
##      B 13 19 25
##      C 15 21 27
##
## , , myDim3 = Z
##
##       myDim2
## myDim1  G  H  I
##      A 19 22 25
##      B 20 23 26
##      C 21 24 27

## Example 3. Transpose one slice from an array.
## http://stackoverflow.com/questions/13811133/apply-a-function-to-each-layer-of-a-3d-array-returning-an-array
A <- array (1:27, c(3,3,3))
A

## , , 1
##
##      [,1] [,2] [,3]
## [1,]    1    4    7
## [2,]    2    5    8
## [3,]    3    6    9
##
## , , 2
##
##      [,1] [,2] [,3]
## [1,]   10   13   16
## [2,]   11   14   17
## [3,]   12   15   18
##
## , , 3
##
##      [,1] [,2] [,3]
## [1,]   19   22   25
## [2,]   20   23   26
## [3,]   21   24   27

## similar idea to previous example. Go through first and second dimensions,
##
apply(A, 1:2, t)

## , , 1
##
##      [,1] [,2] [,3]
## [1,]    1    2    3
## [2,]   10   11   12
## [3,]   19   20   21
##
## , , 2
##
##      [,1] [,2] [,3]
## [1,]    4    5    6
## [2,]   13   14   15
## [3,]   22   23   24
##
## , , 3
##
##      [,1] [,2] [,3]
## [1,]    7    8    9
## [2,]   16   17   18
## [3,]   25   26   27

## The plyr package helps here too. I've seen it do amazing things to
## rearrange repeated observation data from wide to long.

library(plyr)

aaply(A,3,t)

## , ,  = 1
##
##
## X1   1  2  3
##   1  1  4  7
##   2 10 13 16
##   3 19 22 25
##
## , ,  = 2
##
##
## X1   1  2  3
##   1  2  5  8
##   2 11 14 17
##   3 20 23 26
##
## , ,  = 3
##
##
## X1   1  2  3
##   1  3  6  9
##   2 12 15 18
##   3 21 24 27