Where spfun meets bsxfun

matlab

Solution

I think I have found a solution that resolves most of my concerns above. Sometimes when I have the `bsxfun` hammer in my hand, the whole world starts to look like a nail. I forget that some of the simple multiplications that I do with `bsxfun` could be just as easily (and arguably more readably) solved using matrix multiplication. While I don't think that my solution to this problem is any more readable, it is orders of magnitude more efficient than my last solution.

% 'This happens once, outside the loop, since the size'
% 'and sparsity structure of A dont change in the loop'
[M N] = size(A);
[row col] = find(A);

%% 'Inside iterative loop'

    % 'Get the product of non-zero row elements'
    row_prod = accumarray(row,nonzeros(A),[],@prod);

    % 'Use row products to compute 'leave-one-out' row products'
    B = spdiags(row_prod,0,M,M)*spfun(@(x) 1./x, A);

I would still be interested in hearing other suggestions if this can be improved.

Problem

Is there way to perform element-wise operations like `bsxfun` does, but only operate on the non-zero elements of a sparse matrix? In particular, for each non-zero element in matrix `A` at position `(i,j)` I want to find the product of all of the non-zero elements in the `i`-th row, except for element `(i,j)`. For example if the `i-th` row looks like this: ``` 0 5 3 0 0 4 0 0 ``` The result should look like this: ``` 0 12 20 0 0 15 0 0 ``` The most obvious solution seem to be to take the product of the non-zero elements along each row, and then divide each element back out of the row product. So, in the above example, the row product is `5 x 3 x 4 = 60`, and then I just divide out `5` `3` and `4` at their respective locations. Given a sparse matrix `A`, this is my best solution so far: ``` [M N] = size(A); [row col vals] = find(A); row_prod = accumarray(row,vals,[],@prod); B = bsxfun(@ldivide, A, row_prod); B = sparse(row,col,B(sub2ind([M N],row,col)),M,N); ``` The first three lines achieve what I want: a column vector representing the product of the non-zero elements of each row. However, there are quite a few problems with the last two lines. - `bsxfun` is going to return a non-sparse matrix the size of `A` - It is going to waste a LOT of cycles dividing by zero unnecessarily. - The result is a matrix comprised mostly of `Inf` or `-Inf` where I really just want zeros. - I can't mask out the `Inf`s because Matlab defines zero times infinity to be `NaN`. Do I just need to bite the bullet and write a for-loop for this? Or is there another way to approach it?

Original source