4

I have this code currently

s1=seq(0,10,length.out=3)
s2=seq(0,10,length.out=3)
d=outer(s1,s2,`-`)
I=outer(s1,s2,`==`)

However I only want the upper triangular of d and I so I am currently doing

d=d[upper.tri(d,diag=T)]
I=I[upper.tri(I,diag=T)]

Is there a way to make this code faster through incorporating the upper triangular matrix in the outer function? Note that this is a reproducible example, however my code is much bigger and I need to decrease the running time as much as possible.

serv-inc
  • 35,772
  • 9
  • 166
  • 188
raK1
  • 519
  • 1
  • 4
  • 13
  • what is `a` and `b` ? – akrun Mar 22 '16 at 15:32
  • @akrun both s1 and s2 are the same and s=b=s1=s2 – raK1 Mar 22 '16 at 15:32
  • @akrun can you please explain how to use combn in this situation ? – raK1 Mar 22 '16 at 15:35
  • Are these vectors of the same length (in your original data)? Perhaps `combn(seq_along(s1), 2, FUN= function(i) s1[i[1]]- s2[i[2]])` – akrun Mar 22 '16 at 15:36
  • @akrun Thank you , your approach cannot consider the diagonal elements , is there a way to incorporate that into the combn function – raK1 Mar 22 '16 at 15:49
  • Yes that is true, the diagonal elements are always going to be 0 for the first case. – akrun Mar 22 '16 at 15:50
  • If order is not that important, `c(rbind(rep(0,length(s1)),combn(seq_along(s1), 2, FUN= function(i) s1[i[1]]- s2[i[2]])))` – akrun Mar 22 '16 at 15:55

1 Answers1

3

outer internally uses rep to stretch its arguments:

    Y <- rep(Y, rep.int(length(X), length(Y)))
    if (length(X)) 
        X <- rep(X, times = ceiling(length(Y)/length(X)))

You can do the same thing using subsetting, but instead generate only the triangular indices. sequence is a useful helper function in this case:

> i <- sequence(1:3)
> j <- rep(1:3, 1:3)
> d = s1[i] - s2[j]
> I = s1[i] == s2[j]
> d
[1]   0  -5   0 -10  -5   0
> I
[1]  TRUE FALSE  TRUE FALSE FALSE  TRUE
Neal Fultz
  • 9,282
  • 1
  • 39
  • 60
  • Hello @Neal, I think this is slower than `outer` if the lengths of `i` and `j` are large. Do you have any suggestions if that is the case? – Shanks Apr 28 '19 at 14:51
  • I would generally suggest using Rcpp for maximum performance. I would be surprised if `outer` were faster in this case because it does twice as much work, compared to calculating the upper triangle directly, but I suppose it's possible. – Neal Fultz Apr 29 '19 at 16:53