# Generate only unique combinations when input contains duplicates

I have a list with repeated elements, such as

list = {a, a, b, c, c, c}

and I'd like a list of the unique ways to choose 3 elements from it:

{{a, a, b}, {a, a, c}, {a, b, c}, {a, c, c}, {b, c, c}, {c, c, c}}


Alas, "unique" means two different things at once in that sentence, and I can't figure out how to achieve both types of uniqueness simultaneously.

I could use Permutations, whose documentation indicates regarding the input that

Repeated elements are treated as identical.

But I will have many results that differ only by rearrangement, and I do not care about order:

Permutations[list, {3}]

{{a, a, b}, {a, a, c}, {a, b, a}, {a, b, c}, {a, c, a}, {a, c, b}, {a, c, c}, {b, a, a}, {b, a, c}, {b, c, a}, {b, c, c}, {c, a, a}, {c, a, b}, {c, a, c}, {c, b, a}, {c, b, c}, {c, c, a}, {c, c, b}, {c, c, c}}


To eliminate the rearrangements, I could try using Subsets instead, but per its documentation,

Different occurrences of the same element are treated as distinct.

As a result I get many duplicate results that I don't want due to the repeated elements of list:

Subsets[list, {3}]

{{a, a, b}, {a, a, c}, {a, a, c}, {a, a, c}, {a, b, c}, {a, b, c}, {a, b, c}, {a, c, c}, {a, c, c}, {a, c, c}, {a, b, c}, {a, b, c}, {a, b, c}, {a, c, c}, {a, c, c}, {a, c, c}, {b, c, c}, {b, c, c}, {b, c, c}, {c, c, c}}


[Frustrated aside: I can't begin to imagine why Mathematica's permutations-generating function treats repeated list items differently than its combinations-generating function.]

I could eliminate the duplicates from either result, but either way, that still requires calculating the full list of nonunique results as an intermediate step, which I expect to be many orders of magnitude longer than the unique results.

Is it possible to get the result I'm after without having to cull a humongously longer list first to get there? The full problem I am working toward would be a list of 100 elements, ~25 unique elements with multiplicities ranging between 1 and 12, and desired subsets of 7 elements. (100 choose 7) is 16 billion, hence my interest in avoiding computing the full nonunique subset list.

• Can you give ballpark estimates of 1) total length of the list, 2) number of distinct elements, 3) typical number of times an element is repeated; for your intended application? – Marius Ladegård Meyer Apr 9 '20 at 7:32
• The full problem I am working toward would be a list of 100 elements, ~25 unique elements with multiplicities ranging between 1 and 12, and desired subsets of 7 elements. (100 choose 7) is 16 billion, hence my interest in avoiding computing the full nonunique subset list. There are also several sub-problems that interest me if I can't get memory usage down enough to solve the full problem in one blow. – thecommexokid Apr 9 '20 at 15:12
• Why not use DeleteDuplicates on your input and called Subsets on that? – b3m2a1 Apr 9 '20 at 16:47
• Thanks for the specifications, I've edited them into the OP :) My solution runs subsets of 7 elements in about two minutes on a list created with Sort@RandomChoice[Range[1,25], 100] whose largest Tally was 9. There were approximately 1.8 million subsets generated. – Marius Ladegård Meyer Apr 9 '20 at 16:52
• Not to toot my own horn, but I've added two improvements to my post that might be of interest :) Importantly, many approaches in the answers here, including my own and the accepted answer, spend a significant amount of their time converting from one representation (how many elements are selected of each kind) to the representation in the OP (list out the selected elements). This conversion is costly simply because there are so many subsets. OP should consider if this is necessary :) – Marius Ladegård Meyer Apr 12 '20 at 11:23

This is a quick-and-dirty, I'm sure there's more in it, busy with other things but will revisit when time permits.

fn[list_, len_] := Module[{t = Tally[list], u = Union[list]},
Flatten[ConstantArray @@@ Transpose[{u, #}]] & /@
Select[Join @@
Permutations /@
IntegerPartitions[len, {Length@Union@list}, Range[0, len]],
And @@ GreaterEqual @@@ Transpose[{t[[All, 2]], #}] &]];


Tested against what appear to be the current speed champions shows improvement:

list = Sort@Mod[Range[300], 17];
AbsoluteTiming[ans1 = subs[list, 8];] // First
AbsoluteTiming[ans2 = pickDistinct[list, 8];] // First
AbsoluteTiming[ans3 = fn[list, 8];] // First

Sort[Sort /@ ans1] == Sort[Sort /@ ans2] == Sort[Sort /@ ans3]


26.6541

29.7349

18.3133

True

• If you want speed, ciao usually provides! +1 – Marius Ladegård Meyer Apr 9 '20 at 16:53
• Just to make sure I'm understanding the approach: You compute all possible unique multisets of the distinct elements, then reject the ones that exceed the actual multiplicities? – thecommexokid Apr 9 '20 at 18:20
• @thecommexokid - basically. – ciao Apr 9 '20 at 23:59
• Basing a solution around IntegerPartitions was absolutely the way to go for me, so this one gets the Accept for my personal use-case. – thecommexokid Apr 10 '20 at 1:17

My take based on recursion. This approach will be good when the total number of elements is large compared to the number of distinct elements. The dependence on how many elements are picked (length of the subsets) is like that of Pascal's triangle, just like for Subsets.

Wrapper:

pickDistinct[list_, num_] := Block[
{dist, count, elemsLeft, n, picked, result},
{dist, count} = Transpose[Tally[list]];
dist = Reverse[dist[[Ordering[count]]]];
count = Reverse[Sort[count]];
n = Length[dist];
picked = ConstantArray[0, n];
elemsLeft = Table[Total@Drop[count, i], {i, 0, n - 1}];

result = Reap[pickDistinctRec[1, num]][[2, 1]];
Table[
Join @@ Table[ConstantArray[dist[[j]], result[[i, j]]], {j, 1, n}]
, {i, Length[result]}
]
]


Recursion:

pickDistinctRec[pos_, leftToPick_] :=

If[pos == n, picked[[pos]] = leftToPick; Sow[picked],
Do[
picked[[pos]] = m;
pickDistinctRec[pos + 1, leftToPick - m]
, {m, Min[leftToPick, count[[pos]]],
Max[0, leftToPick - elemsLeft[[pos + 1]]], -1}
]
]


For an example with 5 different symbols each occuring 10 times:

list = Sort@Mod[Range[50], 5];
AbsoluteTiming[ans1 = pickDistinct[list, 5];]
AbsoluteTiming[ans2 = Union@Subsets[list, {5}];]


{0.002985, Null}

{1.553072, Null}

so a factor 500 faster. Unfortunately I have 10.0.2.0 so I cannot compare to kglr's algorithm.

UPDATE AND IMPROVEMENT:

The following two tweaks give speedups relative to the original. The original is kept since it has been used for comparison in other posts.

Method 1 (minor):

Do the conversion to OP's subset format in-place instead of looping at the end. Also bail out from recursion as soon as there are no elements left to pick; this is important for the case of small subsets from list of many distinct symbols.

pickDistinct1[list_, num_] := Block[
{dist, count, elemsLeft, n, picked, result},
{dist, count} = Transpose[Tally[list]];
dist = Reverse[dist[[Ordering[count]]]];
count = Reverse[Sort[count]];
n = Length[dist];
picked = ConstantArray[0, n];
elemsLeft = Table[Total@Drop[count, i], {i, 0, n - 1}];

Reap[pickDistinctRec1[1, num]][[2, 1]]
]

pickDistinctRec1[pos_, leftToPick_] :=

If[pos == n, picked[[pos]] = leftToPick;
Sow[Join @@ ConstantArray @@@ Transpose[{dist, picked}]],
Do[
picked[[pos]] = m;
pickDistinctRec1[pos + 1, leftToPick - m]
, {m, Min[leftToPick, count[[pos]]],
Max[0, leftToPick - elemsLeft[[pos + 1]]], -1}
]
]

pickDistinctRec1[pos_, 0] :=

Sow[Join @@
ConstantArray @@@
Transpose[{dist, PadRight[Take[picked, pos - 1], n]}]]


Method 2 ((more) major):

Don't use OP's subset format, but instead return the number of times each symbol has been picked. The function now returns two things: a list of distinct symbols, sorted from largest to smalles frequency, and a list of all the subsets in the format above. NOTE: Whether this is better or not depends a whole lot on how the subsets will be processed after, but that's for the user to decide.

pickDistinct2[list_, num_] :=
Block[{dist, count, elemsLeft, n, picked, result}, {dist, count} =
Transpose[Tally[list]];
dist = Reverse[dist[[Ordering[count]]]];
count = Reverse[Sort[count]];
n = Length[dist];
picked = ConstantArray[0, n];
elemsLeft = Table[Total@Drop[count, i], {i, 0, n - 1}];
{dist, Reap[pickDistinctRec2[1, num]][[2, 1]]}
]

pickDistinctRec2[pos_, leftToPick_] :=
If[pos == n, picked[[pos]] = leftToPick;
Sow[picked],
Do[picked[[pos]] = m;
pickDistinctRec2[pos + 1, leftToPick - m], {m,
Min[leftToPick, count[[pos]]],
Max[0, leftToPick - elemsLeft[[pos + 1]]], -1}]]

pickDistinctRec2[pos_, 0] := Sow[PadRight[Take[picked, pos - 1], n]]


Comparison:

list = Sort @ Mod[Range[300], 17];
First @ AbsoluteTiming[ans = pickDistinct[list, 8];]
First @ AbsoluteTiming[ans1 = pickDistinct1[list, 8];]
First @ AbsoluteTiming[ans2 = pickDistinct2[list, 8];]


35.994123

16.761077

7.780696

Check:

ans === ans1


True

Length[ans] == Length[ans2[[2]]]


True

The approach is now comparable to ciao's:

First @ AbsoluteTiming[fnAns = fn[list, 8];]
Sort[Sort /@ fnAns] === Sort[Sort /@ ans]


20.015753

True

and with

fn2[list_, len_] :=
Module[{t = Tally[list], u = Union[list]},
Select[Join @@
Permutations /@
IntegerPartitions[len, {Length@Union@list}, Range[0, len]],
And @@ GreaterEqual @@@ Transpose[{t[[All, 2]], #}] &]
]

First @ AbsoluteTiming[fnAns2 = fn2[list, 8];]
Sort[fnAns2] === Sort[ans2[[2]]]


8.756621

True

Here is a recursive version that is probably based on similar ideas as the one by @Marius, but looks a little simpler (which is subjective of course) and does not use any mutable state:

ClearAll[subs]
subs[list_List, len_] := Map[
List @@ Flatten[#, Infinity, ll]&,
Flatten @ subs[ll[], Counts[list], len]
]

subs[accum_, _, 0] := accum;

subs[accum_, counts_, left_] :=
With[{fst = First @ Normal @ counts[[{1}]]},
With[{elem = First @ fst, countLeft = Last @ fst},
{
(* Add element, update max count for it *)
subs[
ll[accum, elem],
DeleteCases[ReplacePart[counts, Key[elem] -> countLeft - 1], 0],
left -1
],
(* Skip element *)
Replace[
KeyDrop[counts, elem],
{<||> -> Nothing, rc_ :> subs[accum, rc, left]}
]
}
]
]


It uses linked lists to accumulate the individual sublists, and at each step maintains a partially accumulated list, an association with remaining counts for different elements, and the total number of slots left.

Example:

subs[{a, a, b, c, c, c}, 3]

(* {{a, a, b}, {a, a, c}, {a, b, c}, {a, c, c}, {b, c, c}, {c, c, c}} *)


It seems to be reasonably fast on larger lists, although is probably not the absolute fastest code:

list = Sort @ Mod[Range[50], 5];
subs[list, 5] // Length // AbsoluteTiming

(* {0.005467, 126} *)

• Having Leonid compare his method to mine is such an honor to me :) :) Sadly I cannot test it on my machine because 10.0.2.0 does not implement Nothing. Is there a quick workaround I can use? I tried Replace[DeleteCases[KeyDrop[counts, elem], <||>], {rc_ :> subs[accum, rc, left]}] but MMA still complains about "Part {1} of <||>" not existing... – Marius Ladegård Meyer Apr 9 '20 at 11:14
• @MariusLadegårdMeyer You can use <||> -> Sequence[] in place of <||> -> Nothing, has the same effect and should work on earlier versions of MMA. – Leonid Shifrin Apr 9 '20 at 11:29

An approach based on linear integer programming:

We are given a list of the form $$\{a,a,b,c,c,c\}$$, and a desired length of resulting multiplets $$n$$. Generate the list of distinct elements $$q_i$$ and their multiplicities $$m_i$$. Then solve the constrained system of equations $$\sum x_i = n \\ 0 \leq x_i \leq m_i$$ for the $$x_i$$ over the integers. Each of the resulting solutions for the variables $$\{x_i\}$$ will be correspond to a multiplet of the appropriate length, where the value of $$x_i$$ in each solution corresponds to the multiplicity of the element $$q_i$$ in that multiplet.

n = 3;
list = {a, a, b, c, c, c, c};

(* Create list of distinct elements *)
{distelements, counts} = Transpose[Tally[list]]
(* { {a, b, c}, {2, 1, 3} } *)

(* Create list of dummy variables x_i *)
variables = Array[x, {Length[distelements]}];

(* Open up a can of linear programming *)
soln = soln = variables /.
Solve[Join[{Total[variables] == n},
Thread[0 <= variables <= counts]], variables, Integers]
(* {{0, 0, 3}, {0, 1, 2}, {1, 0, 2}, {1, 1, 1}, {2, 0, 1}, {2, 1, 0}} *)

(* Extract the solutions *)
Flatten[Table[ConstantArray[distelements[[i]], #[[i]]], {i, 1,
Length[distelements]}]] & /@ soln
(* {{c, c, c}, {b, c, c}, {a, c, c}, {a, b, c}, {a, a, c}, {a, a, b}} *)


To see how this scales, I ran this code on a set of 100 randomly-chosen letters of the alphabet (26 distinct elements), with n=7. Mathematica took about 20–30 minutes to return a list of ~3 million subsets on my not-that-powerful laptop.

• Due to the wealth of responses, I didn't get a chance to play with this, but it's a delightful approach to the problem. I think it's the closest conceptual analogue to the clunky thing I would have tried to do if I hadn't had any of these answers to guide me. But much smarter. – thecommexokid Apr 10 '20 at 1:25
ClearAll[kMultiSubsets]

kMultiSubsets = Module[{lst = Sort @ #, k = #2, vars = Array[\[FormalX], Length @ #],
lengths = Tally[Sort @ list][[All, -1]]},
Pick[lst, #, 1] & /@ (vars /. Solve[{Total[vars] == k,
And @@ Join[Thread[0 <= vars <= 1],
And[LessEqual @@ #, 0 <= Total@# <= #2] & @@@
Transpose[{TakeList[vars, lengths], lengths}]]}, vars, Integers])] &;


Examples:

list = {a, a, b, c, c, c};

kMultiSubsets[list, 3]

 {{c, c, c}, {b, c, c}, {a, c, c}, {a, b, c}, {a, a, c}, {a, a, b}}

kMultiSubsets[list, 2]

{{c, c}, {b, c}, {a, c}, {a, b}, {a, a}}

kMultiSubsets[list, 4]

{{b, c, c, c}, {a, c, c, c}, {a, b, c, c}, {a, a, c, c}, {a, a, b, c}}


Note: This approach gives the desired result "without having to cull a humongously longer list, ... but it slower than DeleteDuplicates/Union+ Subsets combination.

Subsets[theSet={a, a, b, c, c, c}, {3}] // Union


Alternative,

Take[#, 3] & /@ Table[RotateLeft[theSet, i], {i, Length[theSet]}]


({{a, b, c}, {b, c, c}, {c, c, c}, {c, c, a}, {c, a, a}, {a, a, b}})

Don't forget to restrict count of element to max subset length (there were no need to do this for the given list. You may want to sort inner lists after construction.

• "I could eliminate the duplicates, but that still requires calculating the full list of nonunique results as an intermediate step, which I expect to be many orders of magnitude longer than the unique results. Is it possible to get the result I'm after without having to cull a humongously longer list first to get there?" – thecommexokid Apr 9 '20 at 6:33
• the second approach does not give the full list of subsets. Take, for example, theSet = {a, a, b, c, c, c, d}. – kglr Apr 9 '20 at 8:33
• Yes, I agree. I think in this case one needs in addition to construct some subsets of list of different elements. – Acus Apr 9 '20 at 8:42

{{a, a, b}, {a, a, c}, {a, b, c}, {a, c, c}, {b, c, c}, {c, c, c}}

One way

 list = {a, a, b, c, c, c}
r = Permutations[list, {3}];
r1 = SortBy[#, ToString[#]] & /@ r;
DeleteDuplicates[r1]


## Timings (any one is free to edit the timing and correct them as needed)

Test 1

ClearAll["Global*"];
list = {a, a, b, c, c, c, d, e, f, g, h, m, n, z, k, r};
kMultiSubsets =
Module[{lst = Sort@#, k = #2, vars = Array[\[FormalX], Length@#],
lengths = Tally[Sort@list][[All, -1]]},
Pick[lst, #, 1] & /@ (vars /.
Solve[{Total[vars] == k,
And @@ Join[Thread[0 <= vars <= 1],
And[LessEqual @@ #, 0 <= Total@# <= #2] & @@@
Transpose[{TakeList[vars, lengths], lengths}]]}, vars,
Integers])] &;

pickDistinct[list_, num_] :=
Block[{dist, count, elemsLeft, n, picked, result}, {dist, count} =
Transpose[Tally[list]];
dist = Reverse[dist[[Ordering[count]]]];
count = Reverse[Sort[count]];
n = Length[dist];
picked = ConstantArray[0, n];
elemsLeft = Table[Total@Drop[count, i], {i, 0, n - 1}];
result = Reap[pickDistinctRec[1, num]][[2, 1]];
Table[Join @@
Table[ConstantArray[dist[[j]], result[[i, j]]], {j, 1, n}], {i,
Length[result]}]];
pickDistinctRec[pos_, leftToPick_] :=
If[pos == n, picked[[pos]] = leftToPick; Sow[picked],
Do[picked[[pos]] = m;
pickDistinctRec[pos + 1, leftToPick - m], {m,
Min[leftToPick, count[[pos]]],
Max[0, leftToPick - elemsLeft[[pos + 1]]], -1}]]


Now

RepeatedTiming[
res = SortBy[#, ToString[#]] & /@ Permutations[list, {4}];
DeleteDuplicates[res];]


RepeatedTiming[
kMultiSubsets[list, 4];]


 RepeatedTiming[Union[Subsets[list, {4}]];]


 RepeatedTiming[pickDistinct[list, 4];]


Test 2

list = {a, a, b, c, c, c, d, e, f, g, h, m, n, z, k, r, k, l, j, x, y,
t, w, q, b, b, b, z};

RepeatedTiming[
res = SortBy[#, ToString[#]] & /@ Permutations[list, {5}];
DeleteDuplicates[res];]

(*9.26*)

RepeatedTiming[
kMultiSubsets[list, 5];]

(*4.666*)

RepeatedTiming[Union[Subsets[list, {5}]];]

(*0.0513*)

RepeatedTiming[pickDistinct[list, 5];]

(* 0.916 *)
`
• "I could eliminate the duplicates, but that still requires calculating the full list of nonunique results as an intermediate step, which I expect to be many orders of magnitude longer than the unique results. Is it possible to get the result I'm after without having to cull a humongously longer list first to get there?" – thecommexokid Apr 9 '20 at 6:33
• Nasser, the last method would not generate the full list of subsets b/c it only takes 5-subsets of consecutive elements. – kglr Apr 9 '20 at 8:21
• @kglr I am a little lost between all these methods :) Please feel free to edit the timing part of my post and fix as needed. – Nasser Apr 9 '20 at 8:24
• I also think that the rotation works in case of 3. Otherwise, probably you in addition need to generate some subsets of lists for rotation. – Acus Apr 9 '20 at 8:30