# How to efficiently find multimap k-itemset frequency

I am doing this to learn association rules. The data d is in format {id, value}, for example:

{{0, 29}, {1, 30}, {1, 31}, {1, 32}, {2, 33}, {2, 34}, {2, 35}, {3, 36}, {3, 37}, {3, 38},
{3, 39}, {3, 40}, {3, 41}, {3, 42}, {3, 43}, {3, 44}, {3, 45}, {3, 46}, {4, 38}, {4, 39},
{4, 47}, {4, 48}, {5, 38}, {5, 39}, {5, 48}}


Bigger data example (task I am using it for has x100 amount data)

I want to find number of id's that contain values {x, y} and P(x|y). Something like:

Contains[d_, x_] := Take[Transpose[Select[d, #[[2]] == x &]], 1][[1]]
ContainsBoth[d_, x_, y_] := Length[Intersection[Contains[d, x], Contains[d, y]]]
ProbabilityOfXGivenY[d_, x_, y_] := ContainsBoth[d, x, y]/Length[Contains[d, y]]


I would like to print out a table or export into a file

As my effort, I can write a brute force solution, that works well for input size of 100 instances.

GetRow[d_, {x_, y_}] := {x, y, ContainsBoth[d, x, y],
100 ProbabilityOfXGivenY[d, x, y] // N,
100 ProbabilityOfXGivenY[d, y, x] // N}
Select[GetRow[d, #] & /@
Select[Tuples[DeleteDuplicates[Transpose[d][[2]]], 2],
(#[[1]] < #[[2]]) &], #[[3]] > 1 &];
TableForm[%, TableHeadings -> {Automatic, {"x", "y", "\[Exists]{x,y}", "P(x|y)", "P(y|x)"}}]


edit: I slightly improved my bruteforce solution:

d = Import["https://raw.github.com/gist/4040530/
fb175e2134559351354572a4aa0554d462cae31b/gistfile1.txt", "Package"];
data = Take[d, 100000];

MultiSet[hm_, key_, value_] := hm[key] = value;
ContainsBoth[hm_, x_, y_] := Length[Intersection[MS[x], MS[y]]]
P[hm_, x_, y_] := 100 (ContainsBoth[hm, x, y]/Length[hm[y]])
GetRow[d_, {x_, y_}] := {x, y, ContainsBoth[d, x, y], P[d, x, y] // N, P[d, y, x] // N}

Timing[
Clear[MS];
MultiSet[MS, #[[1, 2]], #[[All, 1]]] & /@ SplitBy[SortBy[data, Last], Last];
parsedData = {#[[1, 1]], #[[All, 2]]} & /@ SplitBy[SortBy[data, First], First];
nrnodes = Length[pairs = DeleteDuplicates[Flatten[Subsets[#, {2}] & /@
parsedData[[All, 2]], 1]]];
result = GetRow[MS, #] & /@ pairs;
nrnodes]


This completes in 41.091 s on my computer.

Export[FileNameJoin[{NotebookDirectory[], "datax.csv"}], Reverse[SortBy[result, #[[3]] &]]]

-
I suspect your question might be, is there a more efficient method ? – image_doctor Nov 8 '12 at 20:01
Well the Select[Tuples[DeleteDuplicates and Take[Transpose[Select[d, #[[2]] == x &]], 1][[1]] are like the dumbest things ever, but I am not sure how to store the data in Mathematica. – Margus Nov 8 '12 at 20:05

(Edited to respond to comments: the algorithm has been changed and evidence added to show how it reproduces the example in the question.)

Why not exploit built-in procedures to do this work? I refer to matrix multiplication (Dot) for the counting and SparseArray for managing its results. These will do fine given that $x$ and $y$ are integers falling within reasonable ranges. To make life easy, I add $1$ to their values (because SparseArray cannot index from $0$) and later subtract $1$--watch for these moves.

Also, it appears that trivial co-occurrences are not to be reported: only table entries where the third column exceeds $1$ should be output.

a = SparseArray[data + 1 -> ConstantArray[1, Length[data]]];
colsums = Total /@ Transpose[a];
x = Union[data[[All, 2]]] + 1;
x = Select[{x, colsums[[x]]}\[Transpose], Last[#] > 1 & ][[All, 1]];
b = a[[All, x]]\[Transpose] . a[[All, x]]


All the needed information is contained in matrix b:

AbsoluteTiming[
results = Flatten[ParallelTable[{x[[i]] - 1, x[[j]] - 1, b[[i, j]],
b[[i, j]] / b[[j, j]], b[[i, j]] / b[[i, i]]},
{i, 1, Length[x] - 1}, {j, i + 1, Length[x]}], 1];]


{66.3937166, Null}

As an example of how this works, let's consider the data in the question:

data = {{0, 29}, {1, 30}, {1, 31}, {1, 32}, {2, 33}, {2, 34}, {2, 35}, {3, 36}, {3, 37}, {3, 38}, {3, 39}, {3, 40}, {3, 41}, {3, 42}, {3, 43}, {3, 44}, {3, 45}, {3, 46}, {4, 38}, {4, 39}, {4, 47}, {4, 48}, {5, 38}, {5, 39}, {5, 48}};


The a and b matrices look like these:

ArrayPlot[a]


ArrayPlot[b]


A tabulation of results shows how the output in the question is correctly reproduced:

TableForm[results, TableHeadings -> {Range[Length[results]], {"x", "y", "Pairs", "x|y", "y|x"}}]


-
Tally does not seem to do any reduction of the data, i.e. Length[data] == Length[Tally[data]]. – Sasha Nov 8 '12 at 23:39
@Sasha Try it--look at the output. Or experiment with a subset of the data. I can't speak to what's in the data, of course: if it happens to be the case that each ordered pair occurs just once, you can eliminate Tally. For the general problem as stated (and illustrated), tallying seems like an essential step. – whuber Nov 8 '12 at 23:41
I see, the intent of Tally was to help get to the SparseArray rules. You could have done SparseArray[Thread[(1+data)->1]] instead. I'm not nitpicking, just trying to understand what you did. – Sasha Nov 8 '12 at 23:47
@Sasha That's fine--I really appreciate the review, because I'm concerned I may have overlooked something crucial (my solution looks too simple!) For an example of why Tally might be needed, consider processing an artificial dataset such as data = RandomInteger[MultivariatePoissonDistribution[30, {20, 40}], 100000]. Note how the output of Tally (in variable t) is needed at the end. – whuber Nov 8 '12 at 23:55
That is efficient code. However it is completely wrong or I miss to understand what you did. Maybe this will help, in case of data = Take[data, 1000]; the result I am looking for looks like: img831.imageshack.us/img831/8020/newsql.png. – Margus Nov 9 '12 at 1:38

data = Import[
"https://raw.github.com/gist/4040530/\
fb175e2134559351354572a4aa0554d462cae31b/gistfile1.txt", "Package"];


Parse the list of {id,value} into the aggregated list {{id1, valuesList1},{id2,valuesList2},...}}.

parsedData = {#[[1, 1]], #[[All, 2]]} & /@
SplitBy[SortBy[data, First], First];


We can now rebuild your function GetRow to compute the count and conditional probabilities in one go for efficiency:

oGetRow[pd_, {x_, y_}] :=
Cases[pd, {id_, vals_} :>
Block[{h1 = MemberQ[vals, x],
h2 = MemberQ[vals, y]}, {Boole[h1], Boole[h2]} /; h1 || h2]],
both},
If[both > 1,
{x, y, both, 100. both/Count[mask[[All, 2]], 1],
{x, y, both, 0, 0}
]
]


Selecting pairs of values, where the count is greater than 0 can be done much more efficiently:

In[14]:= Length[
pairs = DeleteDuplicates[
Flatten[Subsets[#, {2}] & /@ parsedData[[All, 2]], 1]]]

Out[14]= 565322


You original attempt would need to go through approximately 3*10^8 pairs, vast majority of which would need to be discarded.

Now going through the entire data set takes some time, so it can use ParallelMap:

AbsoluteTiming[
res = Cases[
ParallelMap[oGetRow[parsedData, #] &, pairs], {_, _,
Except[0 | 1], _, _}];]


Using 4 cores this takes over 15 minutes.

So I decided to get a subset of pairs to estimate the timing.

In[18]:= AbsoluteTiming[
res = Cases[
ParallelMap[oGetRow[parsedData, #] &, pairs[[;; 1000]]], {_, _,
Except[0 | 1], _, _}];]

Out[18]= {42.486000, Null}


Thus the projected time on a 4-core machine is almost 7 hours, so within reach of an overnight computation.

In[21]:= 42.486 565/3600.

Out[21]= 6.66794


Using 16-core machine it finished much faster.

In[5]:= AbsoluteTiming[
res = Cases[
ParallelMap[oGetRow[parsedData, #] &, pairs], {_, _,
Except[0 | 1], _, _}];]

Out[5]= {2271.047723, Null}

In[6]:= Length[res]

Out[6]= 76660

In[7]:= Export["~/res.m", res]

Out[7]= ~/res.m


If you are interested in the output, let me know.

-
It takes me a bit to follow but, ˇLength[pairs = DeleteDuplicates[Flatten[Subsets[#, {2}] & /@ parsedData[[All, 2]], 1]]]ˇ in this case is 565322 and not 5357. – Margus Nov 8 '12 at 22:19
My bad. The result of my evaluation was for the subset of data. I will correct that – Sasha Nov 8 '12 at 22:23
@Margus It looks like processing the entire data-set will require some 7 hours of a 4-core machine (mine is few years old laptop). – Sasha Nov 8 '12 at 23:09