[R] A slightly unorthodox matrix product.
Berry, Charles
ccberry @end|ng |rom uc@d@edu
Sat Aug 4 21:31:31 CEST 2018
> On Aug 4, 2018, at 11:43 AM, Jeff Newmiller <jdnewmil using dcn.davis.ca.us> wrote:
>
> Sometimes a good old for loop performs best, even if it doesn't look sexy:
>
>
Fair enough, but a vectorized solution beats them all (see below).
Also,
[SNIP]
> # Charles
> ans1b <- function( a, b )
> {
The lapply you put here was from Eric's solution:
> xxx <- lapply( seq.int( nrow( A ) )
> , function( i ) {
> A[ i, ] %o% B[ i, ]
> }
This is what I had in mind:
ans1b.corrected <- function( a, b ) {
yyy <- sapply( seq.int( nrow( a ) )
, function( i ) a[ i, ] %o% b[ i, ]
, simplify = "array"
)
zzz <- aperm( yyy, c( 3, 1, 2 ) )
zzz
}
On my system it is slower than a for loop but a lot faster than your benchmark showed with the superfluous code from Eric's solution.
For speed, a vectorized solution is faster than a for loop by a factor of 3 on my laptop:
ans0 <- function(A,B){
nca <- ncol(A)
ncb <- ncol(B)
j.index <- rep(1:nca, times = ncb)
k.index <- rep(1:nca, each = ncb)
res <- array(A[, j.index] * B[, k.index], c(nrow(A), nca, ncb))
res
}
> microbenchmark(
+ res0 <- ans0(A, B),
+ res1 <- ans1(A, B),
+ res1b <- ans1b.corrected(A, B),
+ res2 <- ans2(A, B),
+ res3 <- ans3(A, B)
+ )
Unit: microseconds
expr min lq mean median uq max neval cld
res0 <- ans0(A, B) 13.281 18.4960 21.52723 19.9905 23.4750 61.556 100 a
res1 <- ans1(A, B) 353.121 369.8635 409.77788 381.5840 444.3290 701.256 100 e
res1b <- ans1b.corrected(A, B) 82.816 89.4185 101.52321 95.4275 107.1700 217.357 100 c
res2 <- ans2(A, B) 49.674 54.4825 61.78278 58.7540 65.5265 172.625 100 b
res3 <- ans3(A, B) 317.772 342.4220 392.25065 360.4675 436.2125 602.346 100 d
>
FWIW, if there are dimnames on A and B, sapply( row names(A), ..., simplify="array") preserves them without further ado.
Chuck
More information about the R-help
mailing list