1

I am trying to find the count of diagonal 1s in each 3x3 tile e.g.

0 0 1         1 0 0
0 1 0         0 1 0
1 0 0    or   0 0 1

from the below 15x15 matrix.

set.seed(99)
mat <- matrix(sample(c(0,1), 225, prob=c(0.8,0.2), replace=TRUE), nrow=15)
print(mat)

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

I expect the output to be 2 for the above matrix. Is there a way to do this with a for loop and if statements?

Rebecca
  • 45
  • 5

2 Answers2

1

We could use outer(). For this we write two small vectorized functions, that count the elements of the diagonal of a 3x3 slice of our matrix; if the sum is 3 we have a valid diagonal.

For the counterdiagonal we borrow code from this solution.

counterdiag <- function(M) M[(n<-nrow(M))^2-(1:n)*(n-1)]

Now all we need is some coordinates.

m <- n <- mapply(function(i) i:(i+2), 1:13)

And our counting functions.

fun1 <- Vectorize(function(x, y) sum(diag(mat[m[,x], n[,y]])) == 3, SIMPLIFY=FALSE)
fun2 <- Vectorize(function(x, y) sum(counterdiag(mat[m[,x], n[,y]])) == 3, SIMPLIFY=FALSE)

Usage

sum(unlist(outer(1:13, 1:13, fun1)))  # diagonals
# [1] 1

sum(unlist(outer(1:13, 1:13, fun2)))  # counterdiagonals
# [1] 3
Community
  • 1
  • 1
jay.sf
  • 60,139
  • 8
  • 53
  • 110
  • thanks for contributing! If possible, could you explain to me the mapply function you used as well as where the '1:13' comes from? – Rebecca Oct 27 '19 at 18:53
  • @Rebecca Sure. Actually we also could do `sapply(1:13, function(i) i:(i+2))`, `mapply()` becomes interesting when there are more than one argument, try `?mapply`. The `1:13` is the number of columns of the `m`s and `n`s, i.e. the possible combinations of slices in each dimension of `mat`. – jay.sf Oct 27 '19 at 18:58
  • @Rebecca I believe yes. You could wrap this into a function and use `lapply()`, see `?lapply`. This might be an example of how to proceed: https://stackoverflow.com/a/50631083/6574038 – jay.sf Nov 02 '19 at 10:24
0

Here's a nested for loop (using sapply()). Note I did not have the same dataset as you so there's a different seed.

set.seed(123)
mat <- matrix(sample(c(0,1), 225, prob=c(0.8,0.2), replace=TRUE), nrow=15)

n_by_n <- 3L

reg_diag <- diag(n_by_n)
rev_diag <- reg_diag[nrow(reg_diag):1, ]

sum(
  sapply(seq_len(ncol(mat)- n_by_n + 1),
       function(col) {
         sapply(seq_len(nrow(mat) - n_by_n + 1),
                function(row) {
                  tmp <- mat[row:(row + n_by_n - 1), col:(col + n_by_n - 1)]
                  all(tmp == reg_diag) | all(tmp == rev_diag)
                })
       })
)

#[1] 1

If you are only interested in diagonals and do not care about the other values in a submatrix, this splits the matrix by each diagonal and then calculated a rolling sum to see if they sum up to 3:

library(RcppRoll)

set.seed(99)
mat <- matrix(sample(c(0,1), 225, prob=c(0.8,0.2), replace=TRUE), nrow=15)

n_by_n <- 3

diags <- row(mat)- col(mat)
cross_diags <- row(mat) + col(mat)

#could use data.table::frollsum instead of RcppRoll::roll_sumr)
sum(unlist(lapply(split(mat, diags), RcppRoll::roll_sumr, n_by_n), use.names = F) == n_by_n, na.rm = T)
#[1] 1

sum(unlist(lapply(split(mat, cross_diags), RcppRoll::roll_sumr, n_by_n), use.names = F) == n_by_n, na.rm = T)
# [1] 3

A complete base approach would be:

base_rollr <- function(x, roll) {
 #from user @flodel  
    if (length(x) >= roll)  tail(cumsum(x) - cumsum(c(rep(0, roll), head(x, -roll))), -roll + 1)
}

sum(unlist(lapply(split(mat, cross_diags), base_rollr, n_by_n), use.names = F) == n_by_n, na.rm = T)

See also: Get all diagonal vectors from matrix

And: Consecutive/Rolling sums in a vector in R

Cole
  • 11,130
  • 1
  • 9
  • 24
  • Just to clarify @Cole is 'n_by_n' the number of dimensions? When I run this it errors saying the object isn't found so was just wondering. – Rebecca Oct 27 '19 at 18:51
  • Sorry, I skipped a line. Edited but should be 3 – Cole Oct 27 '19 at 18:54
  • @Rebecca, see edit for another approach that only focuses on diagonals instead of looking at the entire 3x3 submatrix. – Cole Oct 27 '19 at 23:55
  • How would I perform the nested for loop (using sapply()) you mentioned first to a list for instance? Would I have to wrap everything into a function for example, `'diag_function'` and then lapply for example, `'lapply(list, diag_function)'`? – Rebecca Oct 29 '19 at 16:56
  • The ```sapply()``` is a double loop. You would wrap everything with ```lapply(mat_list, function (mat) {...})``` – Cole Oct 29 '19 at 19:53
  • How could I get the middle example you mentioned using the RcppRoll package using a list. I have tried wrapping everything with lapply. However, I am getting the error ` Error in roll_sum_impl(x, as.integer(n), as.numeric(weights), as.integer(by), : Not compatible with requested type: [type=list; target=double].` Can you please help @Cole? :) – Rebecca Oct 30 '19 at 11:12
  • I think you need a different question but it'd be ```lapply(list(mat, mat), function (mat) { sum(unlist(...))})```. More broadly the ```RcppRoll``` would not give expected output. Instead of matching identity matrices, it would also evaluate a 3x3 matrix of 1s as matching your criteria. – Cole Oct 30 '19 at 11:18
  • If my approach or jay.sf's solves your question, you should mark a solution. – Cole Oct 31 '19 at 23:07