Independently solving for rows with mrdivide (B/A)?

Question

When I answered this question, my answer started with: "For starters, you can get rid of all the if-statements by storing all the names in a cell."

The "For starters" bit was because I thought that I could add an edit with a vectorized solution. But when I tried to do so I ran into trouble vectorizing the use of mrdivide (b/a).

My question (marked below in the code) is if it is possible to solve b(z,:)/a(z,:) without using a loop. In other words, to solve b/a independently for each row of the matrices.

person = [98   206    35   114;
          60   206    28    52;
         100   210    31   116;
          69   217    26    35;
          88   213    42   100];

person1 = [93 208 34 107];

allNames = {'Cameron'; 'David'; 'Mike'; 'Bill'; 'Joe'};

n = 5;
a = max(person,repmat(person1,n,1));
b = min(person,repmat(person1,n,1));

for z = 1:5
    percent_error = b(z,:)/a(z,:);          %// Here is my question
    if percent_error >= 0.85
        disp(['Match, its ', allNames{z} ,'!'])
    end
end

I am surprised that it works. What does it do mathematically? (I ask this because there is no such operation in mathematics that "vector divide".) — Barney Szabolcs, Dec 11 '12 at 10:18
@BarnabasSzabolcs Yes, I was surprised too, I have never used it like this before. I am not totally sure about what it does, but I think it minimizes `norm(b*x - a)`. — user1884905, Dec 11 '12 at 10:26
I gave it a whirl. :) yes it does look like you're right :O [I fired up a question about that.](http://stackoverflow.com/q/13818288/1031191) — Barney Szabolcs, Dec 11 '12 at 10:30

Eitan T · Accepted Answer · 2012-12-11T13:43:43.867

You can indeed eliminate the loop by vectorizing the operation. The trick is working with diagonal block matrices. Each block is matrix with only one row (each time a different row). After you create such a block matrix for a and for b, you can use mrdivide:

% # Without loop
tic
A = zeros(size(a) * size(a, 1));
B = zeros(size(b) * size(b, 1));
V = ones(size(a, 2), 1) * (1:size(a, 1));
idx = (0:size(A, 1):numel(A) - 1) + (V(:)' - 1) * size(a, 1) + 1;
A(idx) = a';
B(idx) = b';
X = diag(B / A);
percent_error1 = X(1:size(a, 1):end);
toc

% # With loop
tic
percent_error2 = zeros(5, 1);
for z = 1:5
    percent_error2(z) = b(z,:) / a(z,:);
end
toc

The result is:

Elapsed time is 0.000160 seconds.
Elapsed time is 0.000048 seconds.

percent_error1 =
    0.9741
    0.8516
    0.9670
    0.8221
    0.9611

percent_error2 =
    0.9741
    0.8516
    0.9670
    0.8221
    0.9611

Note that this is one of those cases where matrix division of large arrays takes longer than a for loop.

+1 I had to go through it step by step to understand it. Very clever! Also, very good that you did a run-time check. — user1884905, Dec 11 '12 at 13:29
@user1884905 Thanks. I've added a short explanation to clarify what is being done here. — Eitan T, Dec 11 '12 at 13:44
Yes, I accepted your answer, because I don't think there is an "easier" way of doing this. — user1884905, Dec 11 '12 at 14:46

Barney Szabolcs · Answer 2 · 2012-12-11T11:43:12.140

0

I was thinking about this:

person/person1

But this would only give good result when every index in person is bigger than the corresponding index in person1.

edited Dec 11 '12 at 11:43

answered Dec 11 '12 at 10:45

Barney Szabolcs

11,846
12
66
91

I don't think `a / b` produces the same results as the code in the question. – Eitan T Dec 11 '12 at 11:27
Yepp, but why did it work for me then when i did this: [1 2 3 4; 2 3 4 5]/ [1 2 3 4; 1 2 3 4]? – Barney Szabolcs Dec 11 '12 at 11:31
I get a column of zeroes for `percent_error` when I'm using your method. – Eitan T Dec 11 '12 at 11:33
yes, because I did this actually: [1 2 3 4; 2 3 4 5]/[1 2 3 4] :) oops, sorry, corrected. – Barney Szabolcs Dec 11 '12 at 11:43

Independently solving for rows with mrdivide (B/A)?

2 Answers2