5

Is there a way to use the diag() function in a Matrix without using the built-in function or iteration?

   M<-matrix(1:9, ncol=3) # make a matrix 

    q5b<-function(M){ #function

    }

I know that M[1,1], M[2,2], and M[3,3] will give me the same output as diag(M). However, I can't think of a way to do this without a for loop.

My thought process was I should have a condition where row index == column index in the Matrix then print that value. I appreciate any suggestions.

zx8754
  • 52,746
  • 12
  • 114
  • 209
user3018479
  • 165
  • 2
  • 8

2 Answers2

8

Just subset based on another matrix:

> diag(M)
[1] 1 5 9
> M[matrix(rep(sequence(ncol(M)), 2), ncol = 2)]
[1] 1 5 9

The above would run into a problem in a non-square matrix, so we modify it as below.

As your function, one answer for question 5b could be:

q5b <- function(M) { 
  A <- sequence(ncol(M))[sequence(min(nrow(M), ncol(M)))]
  M[cbind(A, A)]
}

Update: Benchmarks are always fun

library(microbenchmark)

fun1 <- function(M) diag(M)
fun2 <- function(M) M[row(M) == col(M)]
fun3 <- function(M) {
  A <- sequence(ncol(M))[sequence(min(nrow(M), ncol(M)))]
  M[cbind(A, A)]
}    

set.seed(1)
M <- matrix(rnorm(1000*1000), ncol = 1000)

microbenchmark(fun1(M), fun2(M), fun3(M), times = 100)
# Unit: microseconds
#     expr       min        lq     median        uq        max neval
#  fun1(M)  4654.825  4747.408  4822.8865  4912.690   5877.866   100
#  fun2(M) 53270.266 54813.606 55059.0695 55749.062 200384.531   100
#  fun3(M)    66.284    82.321   118.8835   129.361    191.155   100
A5C1D2H2I1M1N2O1R2T1
  • 190,393
  • 28
  • 405
  • 485
8

You can use the functions row and col to find the indices where the column number is identical to the row number:

row(M) == col(M)
#       [,1]  [,2]  [,3]
# [1,]  TRUE FALSE FALSE
# [2,] FALSE  TRUE FALSE
# [3,] FALSE FALSE  TRUE

M[row(M) == col(M)]
# [1] 1 5 9
Sven Hohenstein
  • 80,497
  • 17
  • 145
  • 168
  • Nice solution, but any idea why it's so much slower than `diag` or my alternative? – A5C1D2H2I1M1N2O1R2T1 Dec 10 '13 at 11:52
  • 1
    @AnandaMahto Thanks for the tests. I assume my solution is slower because two integer matrices have to be created whose size is identical to the reference matrix. After this all values in these matrices are compared using `==`. Hence, the solution is not very efficient. – Sven Hohenstein Dec 10 '13 at 12:03