Reputation: 1525
How can one extract the "diagonal" from three-dimensional array in R? For a matrix (2D array) one can use the diag(...)
function. In a similar way, given an N x N x M array, a natural operation is to convert it into an N x M matrix by taking the diagonal from each N x N slice and returning it as a matrix.
It's easy to do this using a loop, but that is not idiomatic R and is slow. Another possibility is to use slightly complex indexing (see my own answer to this question) but it is a bit hard to read. What other alternatives are there? Is there a standard R way to do this?
Upvotes: 2
Views: 2742
Reputation: 68
Another approach is subseting the 3 dimensions array with a 2 dimensions matrix:
a <- array(1:100,dim = c(5,5,4))
ref <- cbind(1:5,1:5,rep(1:4,each= 5))
a[ref]
Output is a vector instead of a matrix. On my computer it is more efficient than apply()
and you can also fill the diagonal values.
Upvotes: 0
Reputation: 263411
Although I'm not enamored of the term "3d.diagonal" for this result, it can be achieved with this simple function (up to identity modulo transpose):
arr <- array(1:27,c(3,3,3) )
apply(arr, 3, function(x) x[row(x)==col(x)] )
# returns same value as diag.3d (arr)
[,1] [,2] [,3]
[1,] 1 10 19
[2,] 5 14 23
[3,] 9 18 27
I think a "real diagonal" would be arr[ cbind(1:3,1:3,1:3) ]
Upvotes: 3
Reputation: 6737
Create an array and fill it by some values:
> a=array(0,c(10,10,5))
> for (i in 1:10) for (j in 1:10) for (k in 1:5) a[i,j,k]=100*i+10*j+k-111
Run the apply
function:
> apply(a,3,diag)
[,1] [,2] [,3] [,4] [,5]
[1,] 0 1 2 3 4
[2,] 110 111 112 113 114
[3,] 220 221 222 223 224
[4,] 330 331 332 333 334
[5,] 440 441 442 443 444
[6,] 550 551 552 553 554
[7,] 660 661 662 663 664
[8,] 770 771 772 773 774
[9,] 880 881 882 883 884
[10,] 990 991 992 993 994
Upvotes: 4
Reputation: 49448
Various diagonals:
A = array(1:12, c(2, 2, 3))
apply(A, 1, diag)
# [,1] [,2]
#[1,] 1 2
#[2,] 7 8
apply(A, 2, diag)
# [,1] [,2]
#[1,] 1 3
#[2,] 6 8
apply(A, 3, diag)
# [,1] [,2] [,3]
#[1,] 1 5 9
#[2,] 4 8 12
Upvotes: 3
Reputation: 1525
One possible approach is to use indexing, where the indices are a matrix with three columns. For example:
diag.3d <- function(A) {
# Expect a N x N x M array
stopifnot(length(dim(A)) == 3)
n <- nrow(A)
stopifnot(n == ncol(A))
m <- dim(A)[3]
IXS <- cbind(1:n, 1:n, rep(1:m, each = n))
cn <- colnames(A)
rn <- dimnames(A)[[3]]
matrix(A[IXS], ncol = n, byrow = T, dimnames = list(rn, cn))
}
Although indices (in variable IXS
) seem hard to read.
Upvotes: 1