15
$\begingroup$

I'm running a Monte Carlo simulation. At the end, I want to count how many simulations return a list of same element.

My problem:

SimulationResult = (* n x m matrix of real numbers *)
c = Count[SimulationResult,r_/;Length[Union[r]] == 1]

n is measured in the millions. m can go up to 50. My benchmark shows the Count function is taking roughly 30% of total running time when n = 10M and m = 20.

Is there anyway I can speed up the counting function?

$\endgroup$
  • 2
    $\begingroup$ Try DeleteDuplicates instead of Union $\endgroup$ – Algohi Nov 15 '14 at 22:57
10
$\begingroup$

If they're machine sized reals, then the following should be fast:

n = 10^6;
m = 50;
res = RandomReal[1, {n, m}];  (* SimulationResult *)
res[[RandomSample[Range[n], 1000]]] = ConstantArray[0., m];

cf = Compile[{{x, _Real, 2}},
   Total[1 - Unitize[Max /@ x - Min /@ x]]];

cf[res] // AbsoluteTiming
(*  {0.250830, 1000}  *)

Variations on the body of cf. No real difference in speed, apparently, although they look like they use fewer operations.

Total[UnitStep[Min /@ x - Max /@ x]]  (* most elegant, perhaps *)
Length[x] - Total[Unitize[Max /@ x - Min /@ x]]
$\endgroup$
  • $\begingroup$ Thanks! You take the cheese. Hard to beat native C code. $\endgroup$ – Code Different Nov 15 '14 at 23:54
  • $\begingroup$ @ Michael E2, Zoff Dino: Sorry, but in spite of much of the effort shown in this thread there's really no need for Compile. The duration of Michael E2's elegant "raw" command is (V8.0 Windows) Total[1-Unitize[Max/@res-Min/@res]]//AbsoluteTiming -> {0.3750215,1000} which is almost exactly the same as for the compiled version ({0.3710212, 1000}). $\endgroup$ – Dr. Wolfgang Hintze Nov 20 '14 at 7:53
  • $\begingroup$ @Dr.WolfgangHintze Yes, I saw the same modest improvement on packed arrays. That's because of the vectorization & optimization of such numeric functions on packed arrays does not improve with compilation, as I took Zoff Dino to imply. The compile cf is faster on unpacked numeric arrays, because it seems to automatically convert the input array to a packed array faster than Developer`ToPackedArray. Simon Woods' improvement, which I briefly considered and stupidly didn't even test, is a more significant improvement since it removes one application of Map. $\endgroup$ – Michael E2 Nov 20 '14 at 12:43
14
$\begingroup$
Count[Equal @@@ SimulationResult, True]
Count[SameQ @@@ SimulationResult, True]
Count[Unitize[Max /@ SimulationResult - Min /@ SimulationResult], 0] (*credit: MichaelE2 *)
Count[Length /@ DeleteDuplicates /@ SimulationResult, 1]
Count[SimulationResult, _?(Equal @@ # &)]
Total[Boole[Equal @@@ SimulationResult]]

Timings using @MichaelE2's setup:

Functions:

ClearAll[a, r, c0, c1, c2, c3, c4, c5, c6, c7, c8, tf, cfr, cfi, cf2r,
        cf2i, jcfr, jcfi, functions, names];
functions = {c0, c1, c2, c3, c4, c5, c6, c7, c8, tf, cfr, cf2r, jcfr, cfi, cf2i, jcfi};
names = ToString /@ functions;

c0 = Count[#, r_ /; Length[Union[r]] == 1] &;
c1 = Count[#, r_ /; Length[DeleteDuplicates[r]] == 1] &;
c2 = Count[#, r_ /; Max[r] == Min[r]] &;
c3 = Count[#, {(a_) ..}] &;
c4 = Count[Equal @@@ #, True] &;
c5 = Count[SameQ @@@ #, True] &;
c6 = Count[Unitize[Max /@ # - Min /@ #], 0] &;
c7 = Count[#, _?(Equal @@ # &)] &;
c8 = Count[Length /@ DeleteDuplicates /@ #, 1] &;
tf = Total[Boole[Equal @@@ #]] &;
cfr = Compile[{{x, _Real, 2}}, Total[1 - Unitize[Max /@ x - Min /@ x]]];
cfi = Compile[{{x, _Integer, 2}}, Total[1 - Unitize[Max /@ x - Min /@ x]]];
cf2r = Compile[{{x, _Real, 2}}, 
         Length[x]-Total[Unitize[Block[{y = #}, Max[y] - Min[y]] & /@ x]]];
cf2i = Compile[{{x, _Integer, 2}}, 
         Length[x]-Total[Unitize[Block[{y = #}, Max[y] - Min[y]] & /@ x]]];
jcfr = Compile[{{reals, _Real, 2}}, 
   Block[{nC = 0, nRow, nCol, firstOfRow, row}, {nRow, nCol} = Dimensions[reals];
    Do[row = reals[[i]]; firstOfRow = row[[1]];
     Do[If[firstOfRow != row[[j]], nC++; Break[];], {j, 2, nCol}], {i,nRow}];
    nRow - nC], CompilationTarget -> "C"];
jcfi = Compile[{{reals, _Integer, 2}}, 
   Block[{nC = 0, nRow, nCol, firstOfRow, row}, {nRow, nCol} = Dimensions[reals];
    Do[row = reals[[i]]; firstOfRow = row[[1]];
     Do[If[firstOfRow != row[[j]], nC++; Break[];], {j, 2, nCol}], {i, nRow}];
    nRow - nC], CompilationTarget -> "C"];

NOTE: I don't have a C compiler installed. So, in Jacob's function are compilation target is reverted to WVM thereby eliminating the real advantage of this function.

Real data:

n = 10^6; m = 50;
res = RandomReal[1, {n, m}];
res[[RandomSample[Range[n], 1000]]] = ConstantArray[0., m];
nms = names[[;; -4]];
fns = functions[[;; -4]];
fns2 = fns; fns2[[-1]] = Short[fns2[[-1]], 2/3];
Grid[Thread[{nms, fns2, First@AbsoluteTiming[#[res]] & /@ fns}]]

enter image description here

Integer data:

n = 10^6; m = 50;
res = RandomInteger[100, {n, m}];
res[[RandomSample[Range[n], 1000]]] = ConstantArray[0, m];
nms = Join[names[[;; -7]], names[[-3 ;;]]];
fns = Join[functions[[;; -7]], functions[[-3 ;;]]];
fns2 = fns; fns2[[-1]] = Short[fns2[[-1]], 2/3];
Grid[Thread[{nms, fns2, First@AbsoluteTiming[#[res]] & /@ fns}]]

enter image description here

$\endgroup$
  • 2
    $\begingroup$ For the curious, the reason cf is slower on an integer array is that the numbers of type Integer are first converted to Real. The question specified real numbers, but for integers, cf should be defined Compile[{{x, _Integer, 2}}, ...]. For either type, the uncompiled function Total[1 - Unitize[Max /@ # - Min /@ #]]& is only slightly slower than the compiled function on packed arrays. One could also define a function that would call the appropriate type compiled function depending on the input. $\endgroup$ – Michael E2 Nov 16 '14 at 3:13
  • $\begingroup$ Would you consider including my answer? $\endgroup$ – Jacob Akkerboom Nov 18 '14 at 13:03
  • $\begingroup$ @Jacob, just added your and Simon's answers. $\endgroup$ – kglr Nov 18 '14 at 14:10
  • $\begingroup$ @Jacob, since i don't have a C compiler installed, the timing numbers are actually for default compilation target WVM:( $\endgroup$ – kglr Nov 18 '14 at 14:24
  • $\begingroup$ @kguler Alright, I guess it is fine this way. I edited my own answer with the comparison in the CompilationTarget->"C" case. Thank you for including me :). $\endgroup$ – Jacob Akkerboom Nov 18 '14 at 15:44
13
$\begingroup$

Funnily enough, this is pretty fast for a top-level code:

Count[SimulationResult, {(a_) ..}]
$\endgroup$
6
$\begingroup$

I can shave about 15% off Michael E2's time with slight tweaks to the compiled code:

cf2 = Compile[{{x, _Real, 2}}, 
  Length[x] - Total[Unitize[Block[{y = #}, Max[y] - Min[y]] & /@ x]]];
$\endgroup$
  • $\begingroup$ +1. Nice. I thought I considered trying something like this but stupidly dismissed it. Probably didn't like introducing a new variable. But getting rid of a Map now seems an obvious thing to try. $\endgroup$ – Michael E2 Nov 20 '14 at 13:01
5
$\begingroup$

The code in this answer is about 4x faster than code in other answers, and 3x if you improve other answers by setting CompilationTarget-> "C".

My function

Here is my function

jacobCCf =
 Compile[
  {{reals, _Real, 2}},
  Block[
   {nC = 0, nRow, nCol, firstOfRow, row}
   ,
   {nRow, nCol} = Dimensions[reals];
   Do[
    row = reals[[i]];
    firstOfRow = row[[1]];
    Do[
     If[firstOfRow !=  row[[j]]
      ,
      nC++; Break[];
      ]
     ,
     {j, 2, nCol}
     ]
    ,
    {i, nRow}
    ];
   nRow - nC
   ]
  , CompilationTarget -> "C"
  ]

Setting up the comparison

Definitions for comparison, especially functions by others

jacobWVMCf =
  Compile[
   {{reals, _Real, 2}},
   Block[
    {nC = 0, nRow, nCol, firstOfRow}
    ,
    {nRow, nCol} = Dimensions[reals];
    Do[

     firstOfRow = reals[[i, 1]];
     Do[
      If[firstOfRow !=  reals[[i, j]]
       ,
       nC++; Break[];
       ]
      ,
      {j, 2, nCol}
      ]
     ,
     {i, nRow}
     ];
    nRow - nC
    ]
   ];
cfr = Compile[{{x, _Real, 2}}, 
   Total[1 - Unitize[Max /@ x - Min /@ x]]];
cf2 = Compile[{{x, _Real, 2}}, 
   Length[x] - Total[Unitize[Block[{y = #}, Max[y] - Min[y]] & /@ x]]];
cfrC = Compile[{{x, _Real, 2}}, 
   Total[1 - Unitize[Max /@ x - Min /@ x]], 
   CompilationTarget -> "C"];
cf2C = Compile[{{x, _Real, 2}}, 
   Length[x] - 
    Total[Unitize[Block[{y = #}, Max[y] - Min[y]] & /@ x]], 
   CompilationTarget -> "C"];

Functions for comparing and formatting timings

timing = Function[Null, First@Timing@#, HoldAll];
timingAndName = 
  Function[Null, {timing@#, ToString@Unevaluated@#}, HoldAll];
timingsAndNamesTable = 
  Function[Null, 
   TableForm[timingAndName /@ Unevaluated[{##}], 
    TableHeadings -> {None, {"Timings", "Code executed"}}], HoldAll];

formattedTimingsAndComparison = 
  Function[Null, 
   Block[{timingTable}, timingTable = timingsAndNamesTable@##;
    Column[{StringForm["Comparison for n = ``", n], timingTable, 
      If[SameQ @@ resultNames, "results are equal", 
       Row[{"results are ", Style["not ", Bold], "equal"}]]}, 
     Spacings -> 2]], HoldAll];

compare[] := 
 ReplacePart[
  List @@
   MapIndexed[
    Function[Null, Join[resultNames[[# ;; # &@First@#2]], Hold@#@res],
      HoldAll],
    funcNames
    ]
  ,
  {0 -> formattedTimingsAndComparison, {_, 0} -> Set}
  ]

Actual comparison

Now

n = 10^6; m = 50;
res = RandomReal[1, {n, m}];
res[[RandomSample[Range[n], 1000]]] = ConstantArray[0., m];

resultNames = Hold[jacobN, michaelCN, simonCN, jacobWVMN, michaelN, simonN];
funcNames = Hold[jacobCCf, cfrC, cf2C, jacobWVMCf, cfr, cf2];
compare[]

gives

Comparison for n = 1000000

"Timings"   "Code executed"
--------------------------------
0.083284    jacobN = jacobCCf[res]
0.305532    michaelCN = cfrC[res]
0.224818    simonCN = cf2C[res]
0.302028    jacobWVMN = jacobWVMCf[res]
0.427718    michaelN = cfr[res]
0.314330    simonN = cf2[res]  

results are equal

$\endgroup$
1
$\begingroup$

not sure about this but try it.

Count[SimulationResult, r_ /; Max[r] == Min[r]]
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.