outer() equivalent for non-vector lists in R

Question

I understand how outer() works in R:

> outer(c(1,2,4),c(8,16,32), "*")

     [,1] [,2] [,3]
[1,]    8   16   32
[2,]   16   32   64
[3,]   32   64  128

It basically takes 2 vectors, finds the crossproduct of those vectors, and then applies the function to each pair in the crossproduct.

I don't have two vectors, however. I have two lists of matrices:

M = list();

M[[1]] = matrix(...)
M[[2]] = matrix(...)
M[[3]] = matrix(...)

And I want to do an operation on my list of matricies. I want to do:

outer(M, M, "*")

In this case, I want to take the dot product of each combination of matrices I have.

Actually, I am trying to generate a kernel matrix (and I have written a kernel function), so I want to do:

outer(M, M, kernelFunction)

where kernelFunction calculates a distance between my two matrices.

The problem is that outer() only takes "vector" arguments, rather than "list"s etc. Is there a function that does the equivalent of outer() for non-vector entities?

Alternately, I could use a for-loop to do this:

M = list() # Each element in M is a matrix

for (i in 1:numElements)
{
   for (j in 1:numElements)
   {
      k = kernelFunction(M[[i]], M[[j]])
      kernelMatrix[i,j] = k;
   }
} 

but I am trying to avoid this in favor of an R construct (which might be more efficient). (Yes I know I can modify the for-loop to compute the diagonal matrix and save 50% of the computations. But that's not the code that I'm trying to optimize!)

Is this possible? Any thoughts/suggestions?

Answer 1

The outer function actually DOES work on lists, but the function that you provide gets the two input vectors repeated so that they contain all possible combinations...

As for which is faster, combining outer with vapply is 3x faster than the double for-loop on my machine. If the actual kernel function does "real work", the difference in looping speed is probably not so important.

f1 <- function(a,b, fun) {
  outer(a, b, function(x,y) vapply(seq_along(x), function(i) fun(x[[i]], y[[i]]), numeric(1)))
}

f2 <- function(a,b, fun) {
    kernelMatrix <- matrix(0L, length(a), length(b))
    for (i in seq_along(a))
    {
       for (j in seq_along(b))
       {
          kernelMatrix[i,j] = fun(a[[i]], b[[j]])
       }
    }
    kernelMatrix
}

n <- 300
m <- 2
a <- lapply(1:n, function(x) matrix(runif(m*m),m))
b <- lapply(1:n, function(x) matrix(runif(m*m),m))
kernelFunction <- function(x,y) 0 # dummy, so we only measure the loop overhead

> system.time( r1 <- f1(a,b, kernelFunction) )
   user  system elapsed 
   0.08    0.00    0.07 
> system.time( r2 <- f2(a,b, kernelFunction) )
   user  system elapsed 
   0.23    0.00    0.23 
> identical(r1, r2)
[1] TRUE