Finding a subset of triplets of digits 0-9 such that each digit occurs 40 times in each position in the triplets I am trying to generate a list of digit triplets to specify stimuli in an auditory (speech-in-noise) perception experiment. Each triplet has to have three different digits (i.e., no repetition within triplets). I have recordings of digits 0-9, recorded in each triplet position (i.e., I have 30 digits, 0 in 1st position, 0 in 2nd position, 0 in 3rd position, 1 in 1st position, 1 in 2nd position, etc...). This is for a homogenization task, i.e., a task for figuring out the signal to noise ratio (SNR) for each digit in each position so that the level of individual digits can be adjusted upward or downward to make them equally intelligible at a fixed noise level.

In the homogenization task, stimuli will be presented at 10 SNRs, and I want each digit in each position to be presented four times at each SNR, but in different triplets. Hence, I want 4 x 10 x 10 triplets total. (The four repetitions is somewhat arbitrary, but it's based on previous work, and 400 trials works out nicely with respect to the amount of time listeners will spend doing the task.)

There are 720 total triplets. I use Python, so it's easy to get the full list with:

all_triplets = list(itertools.permutations(np.arange(10),3))


I'm finding it very difficult to get a subset of 400 of the 720 possible triplets such that each digit occurs in each position exactly 40 times. That is, I want 400 triplets such that

np.sum(triplets[:,c]==d) == 40


is True for all c in [0,1,2] and all d in [0,1,2,3,4,5,6,7,8,9], and with triplets a 400 x 3 numpy array containing the subset.

I think that such subsets exist, but I am not sure, since I haven't buckled down to try to work out the math (I may not have the ability to do so even if I do buckle down, and I definitely don't have the time to do so, regardless). I can't figure out how to construct any such subsets from scratch (I wrote some aggressively mediocre code to do this that kind of works sometimes) or to select the right elements from the full set (assuming there are some).

I wrote some less mediocre code to sample randomly from the full subset and score each subset based on how much it deviates from the 40-of-each-digit criterion, and I found a subset with these tallies (rows are triplet position, columns are digits):

[[43, 41, 39, 40, 37, 40, 42, 42, 37, 39],
 [43, 41, 39, 38, 41, 39, 44, 38, 40, 37],
 [43, 38, 40, 38, 38, 41, 42, 44, 39, 37]]


Of course, randomly sampling from approximately 1.92e+213 total subsets of size 400 from the 720 total set of triplets isn't terribly efficient.

So, is there a solution to this problem? A commenter on reddit's intuition is that this problem is NP-complete, but I don't know enough about the relevant math or computer science to have particularly strong intuitions about this possibility. Should I just be happy with the best solution I stumbled on?
You could try using integer linear programming (ILP).  Use the following variables:


We have a variable $x_{ijk}$ for each triplet $ijk$, which is 1 if triplet $ijk$ is selected to be included in your set of 400, or 0 otherwise.
We have a variable $y_{ip}$ to count the number of times that digit $i$ appears in position $p$ among the selected triplets.


Now add constraints to capture all of your requirements:


We want to select 400 triplets: $\sum_{i,j,k} x_{ijk} = 400$.
We want each digit in each position to appear 40 times: $\sum_{i,p} y_{ip} = 40$.
We want the $x$'s and $y$'s to be consistent:

$$\begin{align}
  y_{i1} &amp;= \sum_{j,k} x_{ijk}\\
  y_{j2} &amp;= \sum_{i,k} x_{ijk}\\
  y_{k3} &amp;= \sum_{i,j} x_{ijk}
  \end{align}$$


Finally, we require $0 \le x_{ijk} \le 1$ and $0 \le y_{ip} \le 1$ for all variables.

Now give the resulting system to an integer linear programming solver, and see if it can find a solution.

If it can't find a solution after a reasonable amount of time, your next step could be to find as small a subset as possible (not necessarily size 400, but hopefully not too much larger than 400).  Thus, replace $\sum_{i,p} y_{ip} = 40$ with $\sum_{i,p} y_{ip} \ge 40$, and instead of the constraint $\sum_{i,j,k} x_{ijk} = 400$, ask the ILP solver to minimize the objective function $\sum_{i,j,k} x_{ijk}$.  That will find a way to ensure each digit is presented at least as many times as desired, using as few trials as possible.