R: How to do this matrix operation without loops or more efficient?

Question

I'm trying to make this operation matrices, multiplying the first column with 2, 3 and 4, the first hold value, and then multiply the second column with 3 and 4, keep the value of the third and multiply the third column with 4. I want to do this without using a "for" loop, wanted to use functions like sapply or mapply. Does anyone have an idea how to do it?

Example with one line:

a[1,1]*(a[1,2], a[1,3], a[1,4]) = 2 4 4 4
a[1,1] a[1,2]*(a[1,3], a[1,4]) = 2 4 16 16 #keep a[1,1] a[1,2] 
a[1,1] a[1,2] a[1,3] a[1,3]*(a[1,4]) = 2 4 16 256 # #keep a[1,1] a[1,2] a[1,3] 

Input:

> a<- matrix(2,4,4) # or any else matrix like a<- matrix(c(1,8,10,1,4,1),3,3)
> a
     [,1] [,2] [,3] [,4]
[1,]    2    2    2    2
[2,]    2    2    2    2
[3,]    2    2    2    2
[4,]    2    2    2    2

Output:

> a
     [,1] [,2] [,3] [,4]
[1,]    2    4    16    256
[2,]    2    4    16    256
[3,]    2    4    16    256
[4,]    2    4    16    256

EDIT: LOOP VERSION

a<- matrix(2,4,4); 
ai<-a[,1,drop=F]; 
b<- matrix(numeric(0),nrow(a),ncol(a)-1); 

i<- 1; 

for ( i in 1:(ncol(a)-1)){ 

  a<- a[,1]*a[,-1,drop=F]; 
  b[,i]<- a[,1]; 

}

b<- cbind(ai[,1],b); 

b

Answer 1

If I understand correctly, what you are trying to do is, starting with a matrix A with N columns, perform the following steps:

Step 1. Multiply columns 2 through N of A by column 1 of A. Call the resulting matrix A1.

Step 2. Multiply columns 3 through N of A1 by column 2 of A1. Call the resulting matrix A2.

...

Step (N-1). Multiply column N of A(N-2) by column (N-1) of A(N-2). This is the desired result.

If this is indeed what you are trying to do, you need to either write a double for loop (which you want to avoid, as you say) or come up with some iterative method of performing the above steps.

The double for way would look something like this

DoubleFor <- function(m) {
    res <- m
    for(i in 1:(ncol(res)-1)) {
        for(j in (i+1):ncol(res)) {
            res[, j] <- res[, i] * res[, j]
        }
    }
    res
}

Using R's vectorized operations, you can avoid the inner for loop

SingleFor <- function(m) {
    res <- m
    for(i in 1:(ncol(res)-1)) 
        res[, (i+1):ncol(res)] <- res[, i] * res[, (i+1):ncol(res)]
    res
}

When it comes to iterating a procedure, you may want to define a recursive function, or use Reduce. The recursive function would be something like

RecursiveFun <- function(m, i = 1) {
    if (i == ncol(m)) return(m)
    n <- ncol(m)
    m[, (i+1):n] <- m[, (i+1):n] * m[, i]
    Recall(m, i + 1) # Thanks to @batiste for suggesting using Recall()!
}

while Reduce would use a similar function without the recursion (which is provided by Reduce)

ReduceFun <- function(m) { 
    Reduce(function(i, m) {
            n <- ncol(m)
            m[, (i+1):n] <- m[, (i+1):n] * m[, i]
            m
        }, c((ncol(m)-1):1, list(m)), right = T)
}

These will all produce the same result, e.g. testing on your matrix

a <- matrix(c(1, 8, 10, 1, 4, 1), 3, 3)
DoubleFor(a)
#      [,1] [,2] [,3]
# [1,]    1    1    1
# [2,]    8   32 2048
# [3,]   10   10 1000
all(DoubleFor(a) == SingleFor(a) & SingleFor(a) == RecursiveFun(a) & 
    RecursiveFun(a) == ReduceFun(a))
# [1] TRUE

Just out of curiosity, I did a quick speed comparison, but I don't think any one of the above will be significantly faster than the others for your size of matrices, so I would just go with the one you think is more readable.

a <- matrix(rnorm(1e6),  ncol = 1e3)
system.time(DoubleFor(a))
#    user  system elapsed 
#  22.158   0.012  22.220 
system.time(SingleFor(a))
#    user  system elapsed 
#  27.349   0.004  27.415 
system.time(RecursiveFun(a))
#    user  system elapsed 
#  25.150   1.336  26.534 
system.time(ReduceFun(a))
#    user  system elapsed 
#  26.574   0.004  26.626