# How to do equality check of a large matrix and get the corresponding index position?

Say I have a matrix A = 1000 x 5, and I want to compare it's first element (uppermost top left) to each of the 5 elements in the first column of matrix B of size 5 x 5. Whenever the first time the equality is flagged to say element of A is less than or equal to the element in B, the corresponding row position of B will be printed, and procedure is repeated with the next element down in A again compared with the same column of B. Similarly, elements in second column A will compare with second column of B, third column of A, compared with third column of B. This is repeated for the whole matrix of 5000 elements in A. So essentially the output matrix C will be of size 1000 x 5, with values ranging from 0 to 5.

Or is this easier in VBA or Matlab?

• Could you rephrase your question a bit? Especially your 2nd sentence I find terribly hard to understand. Perhaps you could add an example? And what did you try before you got stuck? – Sjoerd C. de Vries Nov 1 '12 at 20:31
• @Sjoerd Matrix element C11 would be if A11 <= B11 then print B11 index(1), elseif A11 <= B12 then print B12 index(2)...elseif A11 <=B15 print index(5) else print 0. Then next cell C12 would be if A12 <= B11 print B11 index(1) elseif A12 <= B12 print B12 index(2)...do you see the pattern – sebastian c. Nov 1 '12 at 21:09
• @sebastiancheung Have you considered using Position? – dearN Nov 1 '12 at 22:00

I think this does what you want, though I may have transposed the matrices...

For the example I'll use 3 x 6 instead of 5 x 1000

Create some test data:

{n, m} = {3, 6};
a = RandomReal[{0, 1}, {n, m}];
b = Sort /@ RandomReal[{0, 1}, {n, n}];


Visualise:

MatrixForm[a] MatrixForm[b] This is the code:

c = MapThread[Mod[1 + Total[UnitStep[Outer[Subtract, #1, #2]], {2}], n + 1] &, {a, b}];
MatrixForm[c] • Hi @Simon Woods you are a genius, dont know how to thank you! – sebastian c. Nov 1 '12 at 23:14
• Hi @Simon Woods, is this scalable to A=(1000 x 5), B=(5 x 5) so that C=(1000 x 5)? If so what is the best way to do so? – sebastian c. Nov 1 '12 at 23:42
• @sebastiancheung, the code above works with {n,m} = {5,1000}, but I may have misinterpreted the question. In my code a and b must have the same Length (i.e. same number of rows) so you might need to transpose a. Note that I've compared across rows of b, (i.e. element a[[i,j]] is compared with the elements of b[[i]]). So you might need to transpose b too. I've also assumed that the each row of b is in increasing numerical order. If that's not the case then my approach doesn't work at all. – Simon Woods Nov 2 '12 at 11:28
• To summarise: you are probably better off with kguler's solution! – Simon Woods Nov 2 '12 at 11:29
 bb = RandomInteger[{0, 5}, {5, 5}];
aa = RandomInteger[{0, 5}, {10, 5}];
Grid[{{"a", "b"}, {aa // TableForm, bb // TableForm}},
Spacings -> {5, 2}, Dividers -> All] Using Select do define a mapping between two lists:

ClearAll[selPosF];
selPosF[{ax_, bx_}] := (Function[{xx},
First[Select[Range[Length[bx]], (xx <= bx[[#]]) &, 1] /. {} -> {0}]] /@ ax);
cc = selPosF /@ Thread[Transpose /@ {aa, bb}];
Grid[{{"a", "b", "c"}, TableForm /@ {aa, bb, Transpose[cc]}},
Spacings -> {5, 2}, Dividers -> All] Alternatively, using Pick:

ClearAll[pickPosF];
pickPosF[{ax_,  bx_}] := (Function[{xx},
First[Pick[Range[Length[bx]], (xx <= #) & /@ bx] /. {} -> {0}]] /@ ax);


For both functions you can also use the inside MapThread as follows:

cc = Transpose@ MapThread[pickPosF[Sequence[{#1, #2}]] &, Transpose /@ {aa, bb}]

• Hi @kguler, you are a great Mathematica user! – sebastian c. Nov 1 '12 at 23:18
• @sebastian, thank you:) – kglr Nov 2 '12 at 1:16
• Hi @kguler, I ran through the first two sections then ran MatrixForm[aa] and MatrixForm[bb] which were fine as expected and run selPosF[{aa, bb}] again but now it simply returned -> {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}. Not sure what I did when calling selPosF[{aa, bb}]? – sebastian c. Nov 2 '12 at 14:00
• Sorry mistake, forgot to run the rest of the code. – sebastian c. Nov 2 '12 at 14:25

Another interesting old-ish question that popped into the sidebar (Not even sure why no upvotes on it - it is an interesting question!)

My take:

I refactored Kguler's pickPosF, simply changing to

mypickPosF[{ax_, bx_}] := (Function[{xx}, First[Pick[Range[Length[bx]], UnitStep[bx - xx],
1] /. {} -> {0}]] /@ ax)


This improved performance greatly (about triple in the few tests I did), but surprisingly it still lagged behind selPosF. So I used the latter as my standard.

I next refactored Simon Wood's clever use of Outer to comply with the needs of the OP. Much to my surprise (I really expected the technique to be crushingly fast), it was rather slow. Perhaps I'll revisit the idea at a later time.

My first idea was that this is an ideal place for using a piecewise function:

Map[(pw[x_] = Piecewise[MapIndexed[{First@#2, x <= #} &, #[]]];
pw /@ #[]) &, Transpose[Transpose /@ {bb, aa}]]


This was indeed quick, in the few tests I did, about the same as selPosF. I wondered how smart Mathematica might be in culling "dead" parts of the function (that is, parts that never need to be checked due to ordering of the comparisons). Turns out, not too smart. Simple fix: do a PiecewiseExpand on the generated piecewise function, which eliminates the dead parts. This turns out to be surprisingly slow, dropping performance significantly.

Ah, we know things about the data and requirements, let's do it ourselves:

Map[(With[{z = #[]},
pw[x_] = Piecewise[MapThread[{#2, x <= #} &,
Split[FoldList[Max, z]] // {First /@ #,
Most@Accumulate[Prepend[Length /@ #, 1]]} &]]];
pw /@ #[]) &, Transpose[Transpose /@ {bb, aa}]]


Now we're getting somewhere: About a 5X+ improvement on selPosF.

But there's got to be a better way...:

newmex[list_, limsx_] := Module[{rows = Length@list, lims, mo, pss, spans, cab, ca, t},
lims = With[{z = #}, Split[FoldList[Max, z]] // {First /@ #,
Most@Accumulate[Prepend[Length /@ #, 1]]} &] & /@ limsx;
mo = MapThread[Ordering[Join[#1, #2]] &, {Transpose@list, lims[[All, 1]]}];
pss = SplitBy[Position[UnitStep[Subtract[mo, rows + 1]], 1], First];
spans = Span @@@ Transpose[{Most@Prepend[#, 0] + 1, # - 1}] & /@pss[[All, All, 2]];
cab = ConstantArray[0, rows];
Reap[Scan[(t = #; ca = cab;
MapThread[(ca[[t[[-1]][[#1]]]] = #2) &, t[[;; 2]]]; Sow[ca]) &,
Transpose[{spans, lims[[All, 2]], mo}]]][[2, 1]]]


Much faster.

Notes on the above: I've split out parts of the code so readers can see step-by-step what's up. You'd probably want to compact it, though that won't affect performance.

The function takes the comparison array already transposed. So for the OP, simply transpose the comparison array before calling. Doing it this way allows for ragged column comparison lists, kind of useful (row-padding the comparison array with Infinity could also be used, I suppose).

Quick speed comparison to selPosF on a 200 column random integer array, with comparison array columns generated from binning each column of the target by 10, selecting 50 and sorting them. (The latter forces worst-case behavior in my function): 