# Improving selection of elements from large multidimensional array

Frequently in my profession I deal with large arrays of data. Typically these are mixed types including integers, strings, dates, and reals. To recreate a typical data array:

randstring := StringJoin@FromCharacterCode@RandomInteger[{97, 122}, RandomInteger[{5,20}]]
testdata = RandomInteger[{0, 100000}, {1200000, 10}];
testdata[[All, 10]] = Table[randstring, {1200000}];
testdata[[All, 5]] = testdata[[All, 5]] + RandomReal[];

ByteCount[testdata]/(1024.^2)


367.354

So, we'll assume column 1 contains a unique identifier for each list. It is desired to select all lists with a certain unique identifier and perform some kind of calculation on that. As an example we'll calculate the Quartiles of the data in the 5th column.

selquarts[id_] := With[{sel = Select[testdata, #[[1]] == id &]},
If[sel != {}, Quartiles[sel[[All, 5]]], {}]]

id = 5555;
selquarts[id] // Timing


{1.092007, {28123.4, 45390.9, 65606.4}}

So, ~1.1 seconds to select these items. We can speed it up with Cases:

selquarts2[id_] := With[{sel = Cases[testdata, {id, __}]},
If[sel != {}, Quartiles[sel[[All, 5]]], {}]]

selquarts2[id] // Timing


{0.171601, {28123.4, 45390.9, 65606.4}}

Much better, but knowing there are roughly 100,000 unique ids, calculating this will take ~4.75 hours.

Here are the ways I've sped it up so far:

Create a "position lookup" vector:

ids = testdata[[All,1]];
selquarts3[id_] := With[{pos = Flatten@Position[ids, id]}, If[pos != {}, Quartiles[testdata[[pos, 5]]], {}]]

selquarts3[id] // Timing


{0.062400, {28123.4, 45390.9, 65606.4}}

Create a "conditioned dataset", however if your calculation depends on something other than the data in column 5, it becomes less practical.

conditioned = GatherBy[testdata, First];
conditioned = {#[[1, 1]], #[[All, 5]]} & /@ conditioned;
selquarts4[id_] := With[{sel = Cases[conditioned, {id, {__}}][[1, 2]]},
If[sel != {}, Quartiles[sel], {}]]

selquarts4[id] // Timing


{0.015600, {28123.4, 45390.9, 65606.4}}

So, a long worded question to simply ask are there better ways to do this? Maybe using Dataset or Associations? Emphasis is on speed more than memory usage.

You can use Association very fast and very elegantly!

p = GroupBy[testdata, First];
q[id_] := Quartiles[p[id][[All, 5]]];

KeyExistsQ[p, 5555]
(* True *)
q[5555] // AbsoluteTiming
(* {0.000185, {28502.8, 56197.3, 88398.3}} *)

• Holy crap, GroupBy.
– kale
Sep 24, 2014 at 19:43
• I don't have v10. Could you compare the Timing vs "my" GatherBy/Downvalues answer? Sep 24, 2014 at 20:09
• @belisarius When I tested, your solution was about 10% slower (for element selection). I perform 10000 repetition to obtain measurable timings. Sep 24, 2014 at 20:14
• @belisarius @ybeltukov, My timings are very inconsistent but looks like running through a loop 10,000 times, the timings are pretty similar (maybe 10% slower using GatherBy/DownValues.
– kale
Sep 24, 2014 at 20:18
• @kale Of course this solution is easier and more elegant. If it's in the (up/down) 10% range of the other one I've no doubt about which one is better :) Sep 24, 2014 at 20:23

If you're going to reuse the list in V9, a little preprocessing may help:

t1 = GatherBy[testdata, First];
t2 = t1[[All, 1, 1]];
MapIndexed[(pos@#1 = #2[[1]]) &, t2];
selquartsNew[id_] := If[Head[pos[id]] =!= pos, Quartiles[t1[[pos[id]]][[All, 5]]], {}, {}]

selquartsNew[42] // Timing

(* {0., {17875.9, 34250.9, 57568.4}} *)


Timings below comes from version 10.0.1 under Windows. I will use this timing function:

SetAttributes[timeAvg, HoldFirst]

timeAvg[func_] := Do[If[# > 0.3, Return[#/5^i]] & @@ Timing@Do[func, {5^i}], {i, 0, 15}]


## Pick

Although not as fast as pre-grouping methods Pick quite fast when used well.

First transpose your data so columns are easy (and fast) to address:

tdat = Transpose[testdata];


If your example is a typical operation you should keep your data in this format as rows are easier to work with in Mathematica. It also gives you the ability to pack entire rows (columns) of data.

Then:

Pick[#5, #1, 5555] & @@ tdat

{74009.2, 41782.2, 85430.2, 83496.2, 73818.2, 98949.2, 89902.2, 9615.24, 99992.2, 22660.2}

Pick[#5, #1, 5555] & @@ tdat // timeAvg

0.0399363


## Index

A second level of optimization can be had by sorting the data by the first column, then using one of several methods to find the starting and ending indexes of your selection.

One method that does not incur additional memory usage is a binary search. I shall use Leonid's functions from Finding all elements within a certain range in a sorted list. With compilation and a specified starting point for the second search these would be even faster, but even without these optimizations it is still very fast.

First place the data into a sorted format (I assume the data may be kept like this):

sorted = SortBy[testdata, {First}]\[Transpose];


Then a search function (remember to load Leonid's code!):

find[id_] := With[{sel = sorted[[1]]}, bsearchMax[sel, id] ;; bsearchMin[sel, id]]


Example:

find[5555]

67268 ;; 67277


Timing of complete extraction:

sorted[[5, find @ 5555]] // timeAvg

0.000139777


### Index look-up

An extension of the second method is to create a look-up table with the IDs and indexes. This is similar to the GatherBy/GroupBy methods except that it does not require you to restructure your data into groups, merely access them in groups by index. To build the table we can use PositionIndex in 10.0.1:

asc = Span @@@ PositionIndex[First @ sorted][[All, {1, -1}]];  (* association *)


Other versions might you something like:

dsp =
With[{sel = sorted[[1]]},
sel[[First@#]] -> Span @@ # & /@ GatherBy[Range@Length@sel, sel[[#]] &][[All, {1, -1}]]
] // Dispatch;  (* Dispatch table *)


Now the extraction and timing (using the Dispatch table):

sorted[[5, 5555 /. dsp]]

{74009.2, 41782.2, 85430.2, 83496.2, 73818.2, 98949.2, 89902.2, 9615.24, 99992.2, 22660.2}

 sorted[[5, 5555 /. dsp]] // timeAvg

1.39777*10^-6


Using the Association is just a bit faster:

sorted[[5, asc @ 5555]] // timeAvg

8.7860*10^-7


The look-up table (in either format) is about 20MB:

ByteCount /@ {asc, dsp}

{19347576, 19347672}

• Nice. As always.
– kale
Sep 24, 2014 at 21:26