MATLAB: matrix multiplication with exponentiation

Question

Normal matrix multiplication computes $C=A*B$, such that $C_{ij}=\sum_k{A_{ik}*B_{kj}}$

I want to compute D, such that $D_{ij}=\sum_k{e^{A_{ik}*B_{kj}}}$

Basically I want to exponentiate each multipication result before they are summed into a single cell. How do I code this in MATLAB? Do I have to use for loop?

Answer 1

Thanks everybody who tried to help!

A more thorough discussion took place here

One of the answers in the link was already mentioned by Luis, but has the drawback that it blows up memory use when A and B are large. Say A is MxK, and B is KxN, the intermediate matrix created by bsxfun is MxKxN in size. We can keep it to MxN using a for loop, which according to my timing result did not cause performance to change from the one-liner given by Luis. (Details in the link above)

A= rand(500,1000);
B= rand(1000,400);
C4=zeros(500,400);
for j=1:400
   C4(:,j) = sum(exp(bsxfun(@times, A, B(:,j)')),2);
end

Answer 2

No, but you'll need a 3D array (or something equivalent to it). If you put $A_{ik}B_{kj}$ into location $(i + (nm) j, k)$ of $Q$, where size(A) = [n, p] and size(b)= [p, m], then you can compute $e^Q$ and sum over columns.

You can construct $Q$ from some combination of "repmat" and "diag", I'll bet...but the code will sure be opaque. :(

And here's the code to show you how it's done, alas without comments. The idea is to form a pair of vectors, the first consisting of repeated elements of $A$, the second full of repeating elements of $B$, so that their pairwise product contains $A_{ik} * B_{kj}$ for all possible $ij$ pairs, and all $k$. You can then (at the marked point) exponentiate each of these before summing (I used a product with a matrix of 1s, but "sum" would work just as well) and then reshaping to get $C$.

I'm pretty sure that you can restructure this to avoid about half the transposes, if you're so inclined. I've omitted the obvious semicolons on every line so that you can easily see what's going on by looking at the output. Hint: it helps to write out out $AB$ in the form of a matrix of sums of products, like

[ 2*7+3*4     2*0+3*5]
[  ...               ]
[  ...               ]

to guide you in seeing how this all works.

function testprod()
A = [2 3; 4 5; 6 7]
B = [7 0; 4 5]
n = size(A, 1)
k = size(A, 2)
p = size(B, 2)

v1 = repmat(A, 1, p)'
v1 = v1(:)
v2 = repmat(B(:), n, 1)

h = v1 .* v2
hr = reshape(h, k, [])
% INSERT hr = 2.^ hr  HERE!
s = ones(1, k) * hr

C = reshape(s, k, n)' 
Cp = A * B

Answer 3

You can do that easily with bsxfun:

sum(exp(bsxfun(@times, permute(A, [1 3 2]), permute(B, [3 2 1]))), 3)