# Optimising Table wrapping over a Select

Suppose we have the following lists:

L0={{"a"}, {"a", "h"}, {"a", "d", "k", "r", "v"}, {"a", "b", "c",
"k"}, {"a", "b", "c", "s", "u"}, {"a", "b", "f", "t"}, {"a", "b",
"e", "l", "n", "o"}, {"a", "b", "d", "n", "o"}, {"a", "b", "d", "e",
"n", "o"}, {"a", "b", "d", "e", "m", "n", "o", "t"}};
L1={{"a", "b"}, {"a", "c"}, {"a", "d"}, {"a", "e"}, {"a", "f"}};


The aim is to go through L0 for each pair in L1 and count how many times the given pair appears in L0 and categorise them into a polynomial with respect to lengths. I do as follow:

Table[Total[
x^Map[Length, Select[L0, SubsetQ[ToLowerCase[#1], L1[[i]]] &]]], {i,
Length[L1]}]


which gives:

{2 x^4 + 2 x^5 + 2 x^6 + x^8, x^4 + x^5, 2 x^5 + x^6 + x^8,
2 x^6 + x^8, x^4}


So to explain more clearly the first pair is {"a","b"}, and we see that {"a","b"} appears in elements that have length 4 in L0 two times, two times in elements of size 5 and so on. I wonder how can I make this Table operation faster. This is a sample I'm showing here the main dataset has thousands of elements and there are thousands of curves.

• I believe you've made a typo here: "we see that {"a","b"} appears in elements that have length 2 in L0 two times". Did you mean to say in elements that have length 4? Jul 12, 2019 at 13:49
• Oh yes that's a typo :) thanks Jul 12, 2019 at 13:50
• Not a problem :) Jul 12, 2019 at 13:51
• Does L1 always consist of list of same length? Or is the Length/@ L1 not constant? Jul 12, 2019 at 14:00
• @HenrikSchumacher elements of L1 are always have length =2 Jul 12, 2019 at 14:04

Map[Total[x^Cases[L0, p:{OrderlessPatternSequence[## & @@ #, ___]} :> Length[p]]] &]@L1


{2 x^4 + 2 x^5 + 2 x^6 + x^8, x^4 + x^5, 2 x^5 + x^6 + x^8, 2 x^6 + x^8, x^4}

Somewhat faster:

disPatch = Dispatch[{p : {OrderlessPatternSequence[## & @@ #, ___]} :>
x^Length[p], {__} -> Nothing} & /@ L1];
Total[Replace[L0, disPatch, 1], {2}]


{2 x^4 + 2 x^5 + 2 x^6 + x^8, x^4 + x^5, 2 x^5 + x^6 + x^8, 2 x^6 + x^8, x^4}

To start off with, this has a 5-fold speed improvement on my machine. First, one helper function:

findIntersectionsByLength[a_, l_] :=
Map[If[Intersection[a, #] == a, Length[#], 0] &, l];


This takes a set a and searches for all members of l for which a intersect l[[i]] is a. If there is a match, then it returns the length of the match, otherwise it returns 0, for each member of l.

Then we tally these intersections over L0, delete all of the 0 returns, replace them with the polynomial form, and total them, for each element of L1. The Sort is included because Intersection sorts its results as well, so this ensures that it won't fail because {"a","b"} != {"b","a"}.

Table[Total[
DeleteCases[
Tally[findIntersectionsByLength[Sort[l], L0]], {0, _}] /. {e_?NumericQ,
n_?NumericQ} :> n x^e], {l, L1}]


To test this for larger samples, I generated a large sample L0 as follows:

alphabet =
FromCharacterCode[
List /@ ToCharacterCode["abcdefghijklmnopqrstuvwxyz"]];
L0 = Table[
RandomSample[alphabet, RandomInteger[{1, Length[Alphabet]}]], {i,
1, 1000000}];


Using the same L1 as in the question, I get the following timings:

AbsoluteTiming[
res1 = Table[
Total[x^Map[Length,
Select[L0, SubsetQ[ToLowerCase[#1], L1[[i]]] &]]], {i,
Length[L1]}];]


{34.0179, Null}

AbsoluteTiming[
res2 = Table[
Total[DeleteCases[
Tally[findIntersectionsByLength[Sort[l], L0]], {0, _}] /. {e_?
NumericQ, n_?NumericQ} :> n x^e], {l, L1}];]


{6.74027, Null}

res1 === res2


True

A bit ugly, but also works:

Dot[Power[x,First[#1]]& /@ Tally@(Length/@Cases[L0,Flatten@{___,First@#,___,Last@#,___}]),#[[2]]& /@ Tally@(Length/@Cases[L0,Flatten@{___,First@#,___,Last@#,___}])]& /@ L1


And seems to be quite quick. For the large sample @eyorble generated I get:

AbsoluteTiming[Dot[Power[x,First[#1]]& /@ Tally@(Length/@Cases[L0,Flatten@{___,First@#,___,Last@#,___}]),#[[2]]& /@ Tally@(Length/@Cases[L0,Flatten@{___,First@#,___,Last@#,___}])]& /@ L1]


{0.623584, {0, 0, 0, 0, 0}}

And for @eyorble's algorithm:

{5.83138, {0, 0, 0, 0, 0}}

{2 x^4 + 2 x^5 + 2 x^6 + x^8, x^4 + x^5, 2 x^5 + x^6 + x^8, 2 x^6 + x^8, x^4}

, as required.

I'd like to employ SparseArrayand thus, I convert the characters to integers in Range[1,26]. Actually, I start with integers and convert them to characters (just to be able to test against the original implementation).

n = 10000;
m = 100;
SeedRandom[1234];
LL0 = Table[RandomChoice[1 ;; 26, {RandomInteger[{1, 26}]}], {m}];
LL1 = RandomChoice[1 ;; 26, {n, 2}];
L0 = (FromCharacterCode[Partition[#, 1]]) & /@ (LL0 + 96);
L1 = (FromCharacterCode[Partition[#, 1]]) & /@ (LL1 + 96);


OP's implementation

First@AbsoluteTiming[
result1 =
Table[Total[
x^Map[Length,
Select[L0, SubsetQ[ToLowerCase[#1], L1[[i]]] &]]], {i,
Length[L1]}];
]


11.3748

findIntersectionsByLength[a_, l_] :=  Map[If[Intersection[a, #] == a, Length[#], 0] &, l];

result2 =
Table[Total[
DeleteCases[
Tally[findIntersectionsByLength[Sort[l], L0]], {0, _}] /. {e_?
NumericQ, n_?NumericQ} :> n x^e], {l, L1}]; // AbsoluteTiming


2.9406

My implementation:

cf = Compile[{{len, _Integer, 1}, {idx, _Integer, 1}},
If[Length[idx] == 1,
Most[{0}],
Part[len, Most[idx]]
],
CompilationTarget -> "C",
RuntimeAttributes -> {Listable},
Parallelization -> True,
RuntimeOptions -> "Speed"
];

First@AbsoluteTiming[
len0 = Length /@ LL0;
len1 = Length /@ LL1;
With[{
i = Join @@ LL0,
j = Join @@ ((0 LL0 + 1) Range[Length[LL0]])
},
A0 = SparseArray[Transpose[{i, j}] -> 1, {26, Length[LL0]}]
];
With[{
i = Join @@ ((0 LL1 + 1) Range[Length[LL1]]),
j = Join @@ LL1
},
A1 = SparseArray[Transpose[{i, j}] -> 1, {Length[LL1], 26}]
];
B = With[{A = A1.A0},
ArrayFlatten[{
{
Ramp[SparseArray[A - SparseArray[len1 - 1] Unitize[A]]],
SparseArray[ConstantArray[1, {Length[A], 1}]]
}
}]
];
result3 = Total[x^data, {2}];
]


0.190386

Its result does not coincide with OP's but at least with eyorblade's one:

result1 === result2
result1 === result3
result2 === result3


False

False

True

I don't know where the problem is...

I'd like to point out that 2/3 of the computation time is wasted for doing symbolic manipulations:

result3 = Total[x^data, {2}]; // AbsoluteTiming // First


0.127073

The coefficient arrya can be generated much faster:

cg = Compile[{{n, _Integer}, {l, _Integer}},
Table[n, {l}],
CompilationTarget -> "C",
RuntimeAttributes -> {Listable},
Parallelization -> True
];

SetSystemOptions["SparseArrayOptions" -> {"TreatRepeatedEntries" -> Total}];
First@AbsoluteTiming[
coeffarray = With[{
i = Join @@ cg[Range[Length[data]], Length /@ data],
j = Join @@ data
},
SparseArray[Transpose[{i, j}] -> 1, {n, 26}]
];
]


0.022717

Testing the coefficient array:

coeffarray.Table[x^k, {k, 1, 26}] == result3


True

So, if one decides to settle with the coefficient array, the code can be made about 120 times faster for this input.

# Remarks

I added the all-1-column SparseArray[ConstantArray[1, {Length[A], 1}]]} to B so that cf gets never fed an empty list ({}) as second argument; compiled functions don't like that.