4

I have a matrix that should be symmetric according to theory, but might not be observed as symmetric in my data. I would like to force this to be symmetric by using the maximum of the two compared cells.

test_matrix <- matrix(c(0,1,2,1,0,1,1.5,1,0), nrow = 3)
test_matrix
#>     [,1] [,2] [,3]
#>[1,]    0    1  1.5
#>[2,]    1    0  1.0
#>[3,]    2    1  0.0

It's easy enough to do this with a double loop.

for(i in 1:3){
  for(j in 1:3){
    test_matrix[i, j] <- max(test_matrix[i, j], test_matrix[j, i]) 
   }
}
test_matrix
#>      [,1] [,2] [,3]
#> [1,]    0    1    2
#> [2,]    1    0    1
#> [3,]    2    1    0

But my matrix is larger than $3x3$, and R's problems with loops are well-documented. I'm also interested in making my code as clean as possible. In fact, I considered putting this on code golf, but this is a real problem that I think others might be interested in.

I've seen this one as well as this one, but mine is different in that those op's seemed to actually have a symmetric matrix that just needed reordering, and I have a matrix that I need to change to be symmetric.

Community
  • 1
  • 1
gregmacfarlane
  • 2,121
  • 3
  • 24
  • 53
  • 1
    I had seen that before, but it seemed that the answers usually dealt with machine precision or nearly-symmetric matrices. It's possible in my case that cells that should be equal are wildly different. I don't think the answers apply, though the question could be similar. – gregmacfarlane Apr 24 '15 at 21:12
  • I mostly agree, though if you scroll down, it looks like BondedDust's nice answer is similar to my own. I also like it's idea of using `pmean()` (in place of `pmax()`) which might be better, depending on the reason for your matrices asymmetry. ([Link here](http://stackoverflow.com/questions/18165320/creating-a-symmetric-matrix-in-r), since the earlier comment has now been removed.) – Josh O'Brien Apr 24 '15 at 21:17
  • 1
    BTW, with regards to "problems with loops are well-documented": I too have heard a bit of the historic rationale, but in the last few years there has been *considerable* improvement in loop performance and behavior. Do you have recent and reproducible indicators? (I'm not poking, I'm earnestly curious!) – r2evans Apr 24 '15 at 21:38

1 Answers1

9

You could use pmax(), which returns the element-wise maxima of a pair of vectors.

pmax(test_matrix, t(test_matrix))
#      [,1] [,2] [,3]
# [1,]    0    1    2
# [2,]    1    0    1
# [3,]    2    1    0

It'll work with a pair of matrices, as here, because: (1) in R, matrices are 'just' vectors with attached (dimension) attributes; and (2) the code used to implement pmax() is nice enough to reattach the attributes of it's first argument to the value that it returns.

Josh O'Brien
  • 159,210
  • 26
  • 366
  • 455
  • I had never seen `pmax()` before, but now that I have, I see that it solves all kinds of problems that I've made up difficult workarounds for. Nice. – gregmacfarlane Apr 24 '15 at 21:08
  • If the matrix is large (I mean really really large), this method would not be memory efficient. In those cases, I would prefer OP's solution. – Randy Lai Apr 24 '15 at 21:34
  • @gregmacfarlane and JoshOBrien: not stated in the question, but you might want to include a check for a square matrix; a non-square matrix will still do *most* of what you want but without any warning or error. – r2evans Apr 24 '15 at 21:36
  • @RandyLai: perhaps an `Rcpp` solution would work better? Its sugar would not necessarily be better memory-wise (though a little bit), if I were dealing with large enough matrices (and I often do), I'd likely try an iterative compiled variant. Granted, I don't think I can answer "*how big* before it is a problem". – r2evans Apr 24 '15 at 21:40
  • 1
    @r2evans. It depends how you use `RCpp`. If you write a `RCpp` function which is equivalent to `pmax(test_matrix, t(test_matrix))`, there will be no gain at all. Of course, if you write a loop in `RCpp` does the same job as OP's code does, it would do the job better. – Randy Lai Apr 24 '15 at 22:16