On 05/08/18 16:41, Thomas Jagger wrote:
> Date: Sat, 4 Aug 2018 17:52:37 +1200
>
> From: Rolf Turner <r.tur...@auckland.ac.nz
<mailto:r.tur...@auckland.ac.nz>>
> To: r-help <r-help@r-project.org <mailto:r-help@r-project.org>>
> Subject: [R] A slightly unorthodox matrix product.
> Message-ID: <e0073e07-31b3-20df-c20d-6c565c857...@auckland.ac.nz
<mailto:e0073e07-31b3-20df-c20d-6c565c857...@auckland.ac.nz>>
> Content-Type: text/plain; charset="utf-8"; Format="flowed"
>
>
> Can anyone think of a sexy way of forming following "product"?
>
> Given matrices A and B, both with m rows, form a 3 dimensional array C
> such that:
>
> C[i,j,k] = A[i,j]*B[i,k]
>
> I *think* that the following does what I want. (I keep confusing
> myself, so I'm not sure!)
>
> library(abind)
> xxx <- lapply(1:nrow(a),function(i,a,b){a[i,]%o%b[i,]},a=A,b=B)
> do.call(abind,c(xxx,list(along=3)))
>
> Is there a cleverer way?
>
> cheers,
>
> Rolf Turner
>
> --
> Technical Editor ANZJS
> Department of Statistics
> University of Auckland
> Phone: +64-9-373-7599 ext. 88276
Dear Rolf,
Try the following:
B<-matrix(1:12,3,4)
C<-as.vector(A[,rep(seq(ncol(A)),ncol(B))])*as.vector(B[,rep(seq(ncol(B)),each=ncol(A))])
dim(C) <- c(nrow(A),ncol(A),ncol(B))
#test it on column 2 should return true
all(C[,,2]==A*B[,rep(2,ncol(A))])
#on all columns (sapply returns 9 rows with 3 columns all values are TRUE)
all( sapply(seq(ncol(C)),function(i) (C[,,i]==A*B[,rep(i,ncol(A))]) ) )
Note that it creates the final array by taking advantage of the
column-major ordering in R.
Initially, we create a vector by multiplying elementwise the 2 vectors
internally associated with each matrix,
finally, we generate our 3D array by adding the dimensions attribute,
a vector of 3 elements.
This method should be fairly fast since we are using internal R matrix
addressing, and not multiple function calls required by lapply ().
I hope this helps
Neat, and much better than anything I could have thought of. However
microbenchmark() indicates that the Chuck Berry/Jeff Newmiller solution
is about twice as fast. Not that this speed difference is Any Big Deal,
but.
Thanks for taking an interest in this obscure query of mine.
cheers,
Rolf
--
Technical Editor ANZJS
Department of Statistics
University of Auckland
Phone: +64-9-373-7599 ext. 88276
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
