1 Introduction and motivation
The current (Coronavirus Disease 2019) pandemic is rapidly spreading and significantly impacts healthcare systems. Stay-at-home and social distancing orders enforced in many countries are supporting the control of the disease's spread, while causing turmoil in the economic balance and in social structures [BBD + 20] . Rapid detection of cases and contacts is an essential component in controlling the spread of the pandemic. In the US, the current estimation is that at least 500,000 Covid-19 tests will need to be performed daily in order to successfully reopen the economy [LVLM20]
Unfortunately, as humanity attempts to limit the global COVID-19 infection, prophylactic actions are grandly slowed-down by the severe shortages of COVID-19 testing kits [BEF20].
There are currently two types of tests for -Molecular diagnostic tests that detect the presence of SARS-COV-2 nucleic acids in human samples. A positive result of these tests indicates the presence of the virus in the body. -Serological diagnostic tests that identify antibodies (e.g., IgM, IgG) to SARS-COV-2 in clinical specimens [WKLT20] . Serological tests, also known as antibody tests, can be helpful in identifying not only those who are ill, but also those who have been infected, as antibodies are still present in their blood. This identification may be important for several reasons. First, this test can differentiate those who are immune to the virus and those who are still at risk. Secondly, identifying populations who have antibodies can facilitate research on the use of convalescent plasma in the development of a cure for [FA20].
Work done when the author was with ÉNS and Ingenico Laboratories.
As mentioned, both tests are in very short supply. Governments have taken several measures to work around this shortage: the FDA 4 has become more liberal on the approval of tests via the EUA (Emergency Use Authorization) [FA20] ; in the US and the UK there are attempts to boost the number of locally produced test kits to reach a throughput of 100,000 kits per day. Those efforts cannot, however, be followed by many countries and there remains entire swaths of Africa, Asia and Latin America that are under concrete threat.
To optimize the use of available tests, reduce costs and save time, pool testing (or group testing) can be considered. The concept is credited to Dorfman [Dor43] who wished to test US servicemen for syphilis. In pool testing, multiple samples are mixed, and the resulting 'pool' is tested using the same amount of material or equipment that would have been required to test one individual sample. [Ste57, SG59] are important refinements of Durfman's work. However, when at least one sample in the pool is positive, then the pool test fails. This means that (at least) one sample in the pool is positive, but gives no information about which one. The most naive approach is then to re-test individually each sample, which can be costly and time consuming.
The area has seen numerous developments.
[AJS19] provides a recent survey of the topic and [boo93] is the reference book on the topic.
It is important to distinguish between two types of pool tests: Adaptive tests where the tested samples depend on previously tested ones and Nonadaptive tests, where all the tests are planned in advance.
Pool tests are also classified as either probabilistic or combinatorial. In essence, probabilistic models assume a probability distribution and seek to optimize the average number of tests required to test all the patients. By opposition, combinatorial algorithms seek to minimize the worst-case number of tests when the probability distribution governing the infection is unknown.
This paper deals with adaptive probabilistic tests.
Pool testing has already been used to screen large portions of the population in scarcely-infected areas (or as a best-effort measure, when test availability was low). Pool testing has been successfully used to identify viral diseases, such as HIV [NABB19] , ZIKA [BMBM17] , and INFLUENZA [VMW + 12]. In addition, Pool testing has been suggested as a screening method for routine HCV, HBV, and HIV -1 PCR donors for a blood-bank [RWS99] . In light of the recent pandemic and the urging need to test vast number of subjects, the idea of Pool testing is becoming more and more appealing.
It is currently the official testing procedure in Israel [BAKS + 20], Germany, South Korea [Cw20] , and some US 5 [Haa20] and Indian 6 states 7 [Tod20]. Field research focusing on reducing the number of tests [FTO + 20, AWD20, GG20] did not analyse prior information strategies but instead provided simulation (or small sample) results showing the benefits of pool testing. In most cases, the existing literature only uses pooling as a way to screen the infection in an emerging context, not as a precise approach to identify which individuals are infected and which are not.
We also note projects meant to reduce the amount of work required for pool testing: e.g. the Origami Assays [Woo20] Project, a collection of open source pool testing designs for standard 96 well plates. The Origami XL3 design tests 1120 patients in 94 assay wells.
Yelin et al.
[YAST + 20] demonstrated that pool testing can be used effectively to identify one positive SARS-COV-2 result within 32 samples, and possibly within 64 samples if the cycles are amplified, with an estimated false-negative rate of 10%. [Täu20] uses a strategy consisting in running 'cross batches', where the same individuals are tested several times but in different pools, which eventually leads to positive sample identification. The resulting approach ends up using more tests overall (since it tests every individual more than once) than the strategy proposed in this work and does not exploit prior information. Similarly, Sinnott-Armstrong et al. [SAKH20] suggested to identify low-risk individuals (i.e. asymptomatic and mild cases) and to test them as a pool using a matrix-based method, so as to reduce the number of tests required by up to eight-fold, depending on the prevalence.
It is hence plausible to assume that a successful emergency application of the refined pool testing procedures described in this paper would improve the COVID-19 testing capacity significantly.
This paper departs from the above approaches by assuming the availability of extra information: the a priori probability that each given test is negative. In practice, we may either assume that such probabilities are given, estimated from patient trust metrics, or are learned from past COVID-19 tests. We assume in this work that these probabilities are known.
We show that it is possible to find positive samples in an optimal way, i.e., by performing on average the minimum number of tests. This turns out to be faster than blind divide-and-conquer testing in the vast majority of settings.
A concrete consequence of this research is the design of testing procedures that are faster and more cost-effective.
Before introducing models and general formulae, let us provide the intuition behind our algorithms.
Let us begin by considering the very small case of two samples. These can be tested individually or together, in a pool. Individual COVID-19 testing claims a minimum two units of work-check one sample, then check the other. Pool-checking them requires a minimum of one COVID-19 test. If it is highly probable that both samples are negative, then pool testing is interesting: If both samples are indeed negative, we can make a conclusion after one test and halve the COVID-19 test's cost. However, if that fails, we are nearly back to square one: One of these samples (at least) is positive, and we don't know which one.
In this paper, we identify when to check samples individually, and when to pool-check them insteadincluding all possible generalizations when there are more than 2 samples. We assume that the probability of a sample being positive is known to us in advance. The result is a testing 'metaprocedure' that offers the best alternative to sequential and individual testing.
To demonstrate: the testing procedure that always works is to test every sample individually, one after the other: This gives the 'naive procedure', which always performs 2 COVID-19 tests, as illustrated in Figure 1 . In this representation, the numbers in parentheses indicate which samples are being tested at any given point. The leaves indicate which samples are negative (denoted 1) or positive (denoted 0), for instance the leaf 01 indicates that only the second sample is healthy. Note that the order in which each element is tested does not matter: There are thus 2 equivalent naive procedures, namely the one represented in Figure 1 , and the procedure obtained by switching the testing order of (1) and (2).
Alternatively, we can leverage the possibility to test both samples together as the set {1, 2}. In this case, pooling the pair {1, 2} must be the first step: Indeed, testing {1, 2} after any other test would be redundant, and the definition of testing procedures prevents this from happening. If the test on {1, 2} is negative, both samples are negative and the procedure immediately yields the outcome 11. Otherwise, we must identify which of the samples 1 or 2 (or both) is responsible for the test's positiveness. There are thus two possible procedures, illustrated in Figure 2 .
Intuitively, the possibility that this procedure terminates early indicates that, in some situations at least, only one test is performed, and is thus less costly than the naive procedure. However, in some situations up to three tests can be performed, in which case it is more costly than the naive procedure.
(1)
(2) 00 01
(2) 10 11 Figure 1 . The "naive procedure" for n = 2 consists in testing each entity separately and sequentially.
(1,2)
(1)
(2) Concretely, we can compute how many COVID-19 tests are performed on average by each approach, depending on the probability x 1 that the first sample is positive, and x 2 that the second is positive. To each procedure, naive, pool-left, pool-right, we associate the following polynomials representing the expected stopping time:
It is possible to see analytically which of these polynomials evaluates to the smallest value as a function of (x 1 , x 2 ). Looking at Figure 3 , we use these expectations to define zones in [0, 1] 2 where each algorithm is optimal (i.e. the fastest on average). More precisely, the frontier between zones C and B has equation x 1 = x 2 , the frontier between A and B has equation x 2 = (x 1 − 1)/(x 1 − 2), the frontier between A and C has equation x 2 = (2x 1 − 1)/(x 1 − 1), and the three zones meet at x 1 = x 2 = (3 − √ 5)/2, a well-known cutoff value observed as early as 1960 [Ung60] . Having identified the zones, we can write an algorithm which, given x 1 and x 2 , identifies in which zone of Figure 3 (x 1 , x 2 ) lies, and then apply the corresponding optimal testing sequence. In the specific case illustrated above, three algorithms out of three were needed to define the zones; however, for any larger scenario, we will see that only a very small portion of the potential algorithms will be carefully selected.
Our objective is to determine the zones, and the corresponding testing algorithms, for arbitrary n, so as to identify which samples in a set are negative and which are not, while minimizing the expected number of testing operations.
This section will formalize the notion of a testing procedure, and the cost thereof, so that the problem at hand can be mathematically described. We aim at the greatest generality, which leads us to introduce 'and-tests', a special case of which are samples that can be pool tested.
We consider a collection of n samples. Let [n] denote {1, . . . , n}, and Ω = P([n])\{ / 0}, where P is the power set (ie. P(X) is the set of subsets of X).
A test is a function φ : Ω → {0, 1}, that associates a bit to each subset of Ω .
We focus in this work on the following:
Definition 2 (And-Tests). An and-test φ : Ω → {0, 1} is a test satisfying the following property:
In other terms, the result of an and-test on a set is exactly the logical and of the test results on individual members of that set.
Remark 1. Note that 'or-tests', where ∧ is replaced by ∨ in the definition, are exactly dual to our setting. 'xor-tests' can be defined as well but are not investigated here. Although theoretically interesting by their own right, we do not address the situation where both and-tests and or-tests are available, since we know of no concrete application where this is the case.
Elements of Ω can be interpreted as n-bit strings, with the natural interpretation where the i-th bit indicates whether i belongs to the subset. We call selection an element of Ω .
Definition 3 (Outcome). The outcome F φ (T ) of a test φ on T ∈ Ω is the string of individual test results:
Our purpose is to determine the outcome of a given test φ , by minimizing in the expected number of queries to φ . Note that this minimal expectation is trivially upper bounded by n.
Definition 4 (Splitting). Let T ∈ Ω be a selection and φ be a test. Let S be a subset of Ω . The positive part of S with respect to T , denoted S T , is defined as the set
where the operation ∧ is performed element-wise. This splits S into two. Similarly the complement S ⊥ T = S − S T is called the negative part of S with respect to T .
We are interested in algorithms that find F φ . More precisely, we focus our attention on the following:
Definition 5 (Testing procedure). A testing procedure is a binary tree T with labeled nodes and leaves, such that:
1. The leaves of T are in one-to-one correspondence with Ω in string representation; 2. Each node of T which is not a leaf has exactly two children, (S ⊥ , S ), and is labeled (S, T ) where S ⊆ Ω and T ∈ Ω , such that (a) S ⊥ ∩ S = / 0 (b) S ⊥ S = S (c) S ⊥ = S ⊥ T and S = S T . Remark 2. It follows from the definition 5 that a testing procedure is always a finite binary tree, and that no useless calls to φ are performed. Indeed, doing so would result in an empty S for one of the children nodes. Furthermore, the root node has S = Ω . Consider a testing procedure T , defined as above. T describes the following algorithm. At each node (S, T ), perform the test φ on the selection T of samples. If φ (T ) = 0, go to the left child; otherwise go to the right child. Note that at each node of a testing procedure, only one invocation of φ is performed.
The tree is finite and thus this algorithm reaches a leaf S final in a finite number of steps. By design,
Remark 3. From now on, we will fix φ and assume it implicitly.
Remark 4. We represent a testing procedure graphically as follows: Nodes (in black) are labeled with T , whereas leaves (in blue) are labeled with S written as a binary string. This is illustrated in Figure 4 for n = 3.
This representation makes it easy to understand how the algorithm unfolds and what are the outcomes: Starting from the root, each node tells us which entity is tested. If the test is positive, the right branch is taken, otherwise the left branch is taken. Leaves indicate which samples tested positive and which samples tested negative from now on.
Remark 5. The successive steps of a testing procedure can be seen as imposing new logical constraints. These constraints ought to be satisfiable (otherwise one set S is empty in the tree, which cannot happen). The formula at a leaf is maximal in the sense that any additional constraint would make the formula unsatisfiable. This alternative description in terms of satisfiability of Boolean clauses is in fact strictly equivalent to the one that we gave.
In that case, T is understood as a conjunction T [i]=1 t i , S is a proposition formed by a combination of terms t i , connectors ∨ and ∧, and possibly ¬. The root has S = . The left child of a node labeled (T, S) is labeled S ⊥ T = S ∧ (¬T ); while the right child is labeled S T = S ∧ T . At each node and leaf, S must be satisfiable.
To determine how efficient any given testing procedure is, we need to introduce a probability measure, and a metric that counts how many calls to φ are performed.
We consider the discrete probability space (Ω , Pr). The expected value of a random variable X is classically defined as:
Let T a testing procedure, and let S ∈ Ω be one of its leaves. The length T (S) of T over S is the distance on the tree from the root of T to the leaf S. This corresponds to the number of tests required to find S if S is the outcome of φ . The expected length of a testing procedure T is defined naturally as:
It remains to specify the probabilities Pr(ω), i.e. for any given binary string ω, the probability that ω is the outcome.
If the different tests are independent, we can answer this question directly with the following result:
Then, ∀ω ∈ Ω , Pr(ω) can be written as a product of monomials of degree 1 in x 1 , . . . , x n , where
Thus L T is a multivariate polynomial of degree n with integer coefficients.
In fact, or-tests provide inherently independent tests. Therefore we will safely assume that the independence assumption holds.
Example 1. Let n = 5 and ω = 11101, then
Remark 6. L T is uniquely determined as a polynomial by the integer vector of length 2 n defined by all its lengths: (T ) = ( T (0...0), . . . , T (1...1)).
We have now introduced everything necessary to state our goal mathematically. Our objective is to identify the best performing testing procedure T (i.e. having the smallest L T ) in a given situation, i.e. knowing Pr(ω) for all ω ∈ Ω .
We can now explain how to generate all the testing procedures for a given n ≥ 2. One straightforward method is to implement a generation algorithm based on the definition of a testing procedure. Algorithm 1 does so recursively by using a coroutine. The complete list of testing procedures is recovered by calling FindProcedure(Ω , Ω \ { / 0}).
Input: S ∈ Ω , C ∈ Ω . Output: A binary tree.
We implemented this algorithm (in Python, source code available upon request). The result of testing procedure generations for small values of n is summarized in Table 1 . The number of possible testing procedures grows very quickly with n. An informal description of Algorithm 1 is the following. Assuming that you have an unfinished procedure (i.e. nodes at the end of branches are not all leaves). For those nodes S, compute for each T the sets S T and S ⊥ T . If either is empty, abort. Otherwise, create a new (unfinished) procedure, and launch recursively on nodes (not on leaves, which are such that S has size 1).
Algorithm 1 terminates because it only calls itself with strictly smaller arguments. We will discuss this algorithm further after describing some properties of the problem at hand.
Once the optimality zones, and the corresponding testing procedures, have been identified, it is easy to write an algorithm which calls the best testing procedure in every scenario. At first sight, it may seem that nothing is gained from doing so -but as it turns out that only a handful of procedures need to be implemented.
This construction is captured by the following definition:
Definition 6 (Metaprocedure). A metaprocedure M is a collection of pairs (Z i , T i ) such that:
1. Z i ⊆ [0, 1] n , Z i ∩ Z j = / 0 whenever i = j and i Z i = [0, 1] n . 2. T i is a testing procedure and for any testing procedure T ,
A metaprocedure is interpreted as follows: Given x ∈ [0, 1] n find the unique Z i that contains x and run the corresponding testing procedure T i . We extend the notion of expected length accordingly:L M = min i L T i ≤ n.
One way to find the metaprocedure for n, is to enumerate all the testing procedures using Algorithm 1, compute all expected lengths L T from the tree structure, and solve polynomial inequalities.
Surprisingly, a vast majority of the procedures generated are nowhere optimal: This is illustrated in Table 2 . Furthermore, amongst the remaining procedures, there is a high level of symmetry. For instance, in the case n = 3, 8 procedures appear 6 times, 1 a procedure appears 3 times, and 1 procedure appears once. The only difference between the occurrences of these procedures -which explains why we count them several times -is the action of the symmetric group S 6 on the cube (see Appendix B for a complete description).
The metaprocedure for n = 3 cuts the unit cube into 52 zones, which correspond to a highly symmetric and intricate partition, as illustrated in Figures 5, 6, and 7. An STL model was constructed and is available upon request. The large number of suboptimal procedures shows that the generate-then-eliminate approach quickly runs out of steam: Generating all procedures for n = 6 seems out of reach with Algorithm 1 8 . The number of zones, which corresponds to the number of procedures that are optimal in some situation, is on the contrary very reasonable. Proof. By induction on n: There are (n + 1) choices of a root node, P(n) choices for the left child, and P(n) choices for the right child. This gives the recurrence P(n + 1) = (n + 1)P(n) 2 , hence the result.
This number grows rapidly and constitutes a lower bound for the total number of procedures (e.g. for n = 8 we have P(n) > 2 184 ). On the other hand, the naive procedure is the one with maximal multiplicity, which yields a crude upper bound C 2 k P(n) on the number of procedures, where C t is the t-th Catalan number defined by:
The zones can be determined by sampling precisely enough the probability space. Simple arguments about the regularity of polynomials guarantee that this procedure succeeds, when working with infinite numerical precision. In practice, although working with infinite precision is feasible (using rationals), we opted for floating-point numbers, which are faster. The consequence is that sometimes this lack of precision results in incorrect results on the zone borders -however this is easily improved by increasing the precision or checking manually that there is no sub-zone near the borders. 1  1  1  2  4  3  3  312  52  4  36585024  181  5 8.926 · 10 20 ? 6 2.242 · 10 55 ? Table 2 . Procedures and metaprocedures for some values of n. The number of zones for n = 5 and 6 cannot be determined in a reasonable time with the generate-then-eliminate approach.
We now focus on some of the properties exhibited by testing procedures, which allows a better understanding of the problem and interesting optimizations. This in effect can be used to prune early the generation of procedures, and write them in a more compact way by leveraging symmetries. We consider in this section a testing procedure T . Proof. First notice that for all T , T , and
T . We will denote both by S bb T T . We have the following : If there exists k, T 1 , . . . , T k , and β 1 , . . . , β k such that
(Ω )
Remark 7. Proposition 3 indicates that testing procedures are, in general, unbalanced binary trees: The only balanced procedure being the naive one.
Lemma 4. If N is the naive procedure, then for any testing procedure T and for all x 1 , . . . , x n such that
In other terms {∀i ∈ [n], 1 2 ≤ x i ≤ 1} is contained in the naive procedure's optimality zone.
Proof. An immediate corollary of Proposition 3 is that for all i ∈ [n], we have ∂ x i L T (x 1 , . . . , x n ) ≥ 0, where ∂ x i indicates the derivative with respect to the variable x i . Since the native procedure has a constant length, it suffices to show that it is optimal at the point { 1 2 , . . . , 1 2 }. Evaluating the length polynomials at this point gives
Now remember that the naive procedure gives the only perfect tree. It suffices to show that unbalancing this tree in any way results in a longer sum in the equation above. Indeed, to unbalance the tree one needs to: -Remove two bottom-level leaves, turning their root node into a leaf -Turn one bottom-level leaf into a node -Attach two nodes to this newly-created leaf
The total impact on the sum of lengths is +1. Hence the naive algorithm is minimal at { 1 2 , . . . , 1 2 }, and therefore, in the region {∀i ∈ [n], 1 2 ≤ x i ≤ 1}.
Remark 8. This also show that if we assume that the probabilities are supposed uniform (ie. we assume no a priori knowledge) the optimal procedure is the naive one. Therefore we can see that the gain for n = 3 is aprroximately 0, 34 since the optimal procedure in average gives 2, 66. In percentage the gain is 15%. If the probabilites are very low we have a gain of almost 2 which is 3 times faster. As expected it is much more interesting if we think that the samples have a good chance to be negative, which is the case in most real life scenarii.
Lemma 5. If the root has a test of cardinality one, then the same algorithms starting at both sons have same expected stopping time. This applies if the next test is also of cardinal one.
Proof. Without loss of generality we can assume that the test is {1}. We have {0, 1} n {1} = {0b 2 · · · b n |b 2 · · · b n ∈ {0, 1} n−1 } and {0, 1} n ⊥ {1} = {1b 2 · · · b n |b 2 · · · b n ∈ {0, 1} n−1 }. A test T that doesn't test 1 applied on those sets will give the same split for both, and the probability that the test answers yes or no is the same. This is also true for the sets and the tests T such that i is not in T for i in {1, . . . , k}. {0 k b 2 · · · b n |b 2 · · · b n ∈ {0, 1} n−k } and {01 k b 2 · · · b n |b 2 · · · b n ∈ {0, 1} n−k }. A test T such that there exists i in T {1, . . . , k} brings no information for the set of possibilities {01 k b 2 · · · b n |b 2 · · · b n ∈ {0, 1} n−k }, but testing this i is useless for the set {0 k b 2 · · · b n |b 2 · · · b n ∈ {0, 1} n−k }. So we can apply the test T − {1, . . . , k}.
Corollary 1. If the root has a test of cardinal one, then an optimal algorithm can always apply the same test for the right and left child. If this test is also of cardinal one then the property is still true.
This result helps in identifying redundant descriptions of testing procedures, and can be used to narrow down the generation, by skipping over obvious symmetries of the naive procedure (see Figure 8 ). To further accelerate generation we can only keep one representative of each algorithms that have the same expected length for all x i . Lemma 6. If a node labeled T 1 has two children that are both labeled T 2 , then we can interchange T 1 and T 2 without changing the testing procedure's expected length.
Yet another simple observation allows to reduce the set of subsets T at each step:
Lemma 7. Consider a node labeled (T, S ). Assume that there is i ∈ [n] such that, for all S in S , i / ∈ S. Then we can replace T by T ∪ {i}.
Proof. We can easily see that S T = S T ∪{i} and S ⊥ T = S ⊥ T ∪{i} .
Finally we can leverage the fact that the solutions exhibit symmetries, which provides both a compact encoding of testing procedures, and an appreciable reduction in problem size.
Lemma 8. Let σ ∈ S n be a permutation on n elements. If we apply σ to each node and leaf of T , which we can write σ (T ), then L σ (T ) (x 1 , . . . , x n ) = L T (σ (x 1 , . . . , x n )) .
Proof. Note that for any S ∈ Ω and T ∈ Ω \ { / 0} we have σ S T = S σ (T ) and σ S ⊥ T = S ⊥ σ (T ) , where σ operates on each binary string. It follows that for any leaf S, T (S) becomes T (σ (S)) under the action of σ , hence the result.
Lemma 9. Let S be a simplex of the hypercube, T a procedure, E = {σ (T )|σ ∈ S n }, then there exists T 0 in E, such that for all x in S, T 1 in E we have
Moreover we have for all σ in S n , x in σ (S),
Remark 9. The last two propositions allow us to solve the problem on a simplex of the hypercube (of volume 1/n!) such as {p 1 , . . . , p n | 1 ≥ p 1 · · · ≥ p n ≥ 0}.
We examine the following problem: Find the testing procedure T for a given k ≤ n, (p i 1 , . . . , p i k ) ∈ [0, 1] n , and a selection P ⊆ 2 [k] that satisfies:
This can be computed using a dynamic programming technique, by examining the outcome of each possible test that is the root node of the testing procedure T , which gives Algorithm 2.
The same dynamic programming algorithm can also be used to compute the number of testing procedures (including those leading to duplicate polynomials) that exist in a given dimension. It is actually even easier, since there is a huge number of symmetries that can be exploited to count. 9
Definition 7 (Decided point). We say that x is a decided point for S a set of selections if either of the following is true:
In the first case, we will say that x is a positive decided point, and a negative decided point in the second case.
We denote by D + S the set of positive decided points of S , D − S its set of negative decided points, and D S = D + S ∪ D − S its set of decided points.
Algorithm 2: FindOptimal
Output: The optimal testing procedure T at point (p 1 , . . . , p k ) which satisfies S T = S .
1. if k == 0 then return the naive algorithm 2. if |D S | > 0 3.
R ← {{r 1 , . . . , r p } | {u r 1 , . . . , u r p } ∪ D + S } 5.
T ← FindOptimal ( , (p u 1 , . . . , p u ), R) 6.
replace {t 1 , . . . ,t r } by {u t 1 , . . . , u t r } in T 7.
replace { 1 , . . . , r } by {u 1 , . . . , u r } ∪ D + S in T 8. else 9 Indeed, we can apply the algorithm to an even higher dimension than our solution to the given point problem.
if S ⊥ = / 0 or S = / 0 then continue 14.
T ⊥ ← FindOptimal(k, (p 1 , . . . , p k ), S ⊥ ) 15.
T ← FindOptimal(k, (p 1 , . . . , p k 
return the best algorithm in W at point (p 1 , . . . , p n )
Counting the number of algorithms in a given dimension works the same way; the only difference is that there is no need to look at the probabilities, and thus, the resulting Algorithm 3 does fewer recursive calls and is faster. 10
Algorithm 3: CountAlgorithms
Output: The number of testing procedures which satisfy S T = S .
c ← CountAlgorithms(k, (p 1 , . . . , p k ), S ) 13.
c ← c + c c ⊥ 14. return c 7 Enumerating procedures for n = 3 All the procedures for n = 3 that are optimal at some point, up to symmetries, are represented in Figure 9 .
We have introduced the question of optimal pool testing with a priori probabilities, where one is given a set of samples and must determine in the least average number of operations which samples are negative, and which are not. We formalized this problem and pointed out several interesting combinatorial and algebraic properties that speed up the computation of an optimal sequence of operations -which we call a metaprocedure. We determined the exact solution for up to 4 samples.
For larger values, our approach requires too many computation to be tractable, and thus an exact solution is out of reach; however we gave several heuristic algorithms that scale well. We showed that these heuristics are sub-optimal in all cases, but they always do better than standard screening. The existence of a polynomial-time algorithm that finds optimal metaprocedures for large value of n is an open question -although there is probably more hope in finding better heuristics. An alternative would be to modify our generation algorithm to kill branches when the resulting expected lengths are all worse than some already-known procedure.
Once the metaprocedure for a given n is known, which only needs to be computed once, implementation is straightforward and only invokes a handful of (automatically generated) cases. 10 We are not aware of a closed-form formula providing the same values as this algorithm. Finally, in our model we do not consider false positives and false negatives. In other words, tests are assumed to be 100% accurate. Integrating in the model false positive and false negative probabilities is an interesting research challenge.
Besides the performance gain resulting from implementing metaprocedures for sample testing, the very general framework allows for applications in medical and engineering tests.
The approach consisting in generating many candidates, only to select a few, is wasteful. In fact, for large values of n (even from 10), generating all the candidates is beyond reach, despite the optimizations we described.
Instead, one would like to obtain the optimal testing procedure directly. It is a somewhat simpler problem, and we can find the solution by improving on our generation-then-selection algorithm (see Appendix 6). However if we wish to address larger values of n, we must relax the constraints and use the heuristic algorithms described below, which achieve near-optimal results. This would be usefull in real life scenarii for COVID-19 tests since we would like to test hundreds or more samples to have real gain.
We first associate a 'cost' to each outcome S, and set of outcomes S : Intuitively, we give higher gains to subsets T on which testing gives more information. Note that, if a call to φ doesn't give any information (i.e. S T or S ⊥ T is empty), then gain(T, S) = 0. This heuristic provides us with a greedy algorithm that is straightforward to implement. For given values x 1 , . . . , x n we thus obtain a testing procedure T H .
Testing the heuristic. We compared numerically T H to the metaprocedure found by exhaustion in the case n = 3. The comparison consists in sampling points at random, and computing the sample mean of each algorithm's length on this input. The heuristic procedure gives a mean of 2.666, which underperform the optimal procedure (2.661) by only 1%.
Counter-example to optimality. In some cases, the heuristic procedure behaves very differently from the metaprocedure. For instance, for n = 3, x 1 = 0.01, x 2 = 0.17, x 3 = 0.51, the metaprocedure yields a tree which has an expected length of 1.889. The heuristic however produces a tree which has expected length 1.96. Both trees are represented in Figure 10 .
Beyond their different lengths, the main difference between the two procedures of Figure 10 begin at the third node. At that node the set S is the same, namely To understand this difference, first notice that besides 010 and 011, all leaves are associated to a very low probability. The heuristic fails to capture that by choosing T = {1, 3} early, it could later rule out the leaf 010 in one step and 011 in two. There does not seem to be a simple greedy way to detect this early on. (1,2,3)
(1,2)
(1)
(3) The optimal procedure has expected length 1.889, as compared to 1.96 for the heuristic procedure.
Another approach is to use small metaprocedures on subsets of the complete problem. Concretely, given n samples to test, place them at random into k-tuples (from some small value k, e.g. 5). Then apply the k-metaprocedure on these tuples. While sub-optimal, this approach does not yield worst results than the naive procedure.
In cases where it makes sense to assume that all the x i are equal, then we may even recursively use the metaprocedures, i.e. the metaprocedures to be run are themselves places into k-tuples, etc. Using lazy evaluation, only the necessary tests are performed.
B Equivalences and symmetries for n = 3 A procedure can undergo a transformation that leaves its expected length unchanged. Such transformations are called equivalences. On the other hand, Lemma 8 shows that some transformations operate a permutation σ on the variables x i -such transformations are called symmetries.
Equivalences and symmetries are responsible for a large part of the combinatorial explosion observed when generating all procedures. By focusing on procedures up to symmetry, we can thus describe the complete set in a more compact way and attempt a first classification.
In the following representations (Figures 11, 12 , and 13), blue indicates a fixed part, and red indicate a part undergoing some permutation. Double-headed arrows indicate that swapping nodes is possible. The number of symmetries obtained by such an operation is indicated under the curly brace below. Figure 11 . Trees representation with a grouping by one element on the root. For a fixed element, we have 2 2 possible permutations. Since we have 4 patterns, we get 2 2 × 4 possible permutations for one grouping. Hence, we finally have 2 2 × 4 × 3 for all possible groupings by one element.
(1,2)
(1) Figure 12 . Tree representations with a grouping by two elements on the root. For 10 fixed elements, we have 2 possible permutations, for 2 fixed elements, we have 2 possible permutations, and for 2 possible permutations, we have 6 possible permutations. Hence, we finally have 2 × 10 + 4 × 2 + 6 × 2 for all possible groupings by two elements.
(1,2,3)
(1)(2)(3) (1,2,3)
All permutations with g = 2 (3) Figure 13 . Trees representation with a grouping by three elements on the root. For a fixed element at the upper left corner side, we have 2 2 possible permutations. For the upper right corner side, we get 2 2 . We replace the subroot of the fixed trees and get (2 2 + 2 2 ) × 3. We also have the 40 × 3 trees from the grouping of two (g = 2). Hence, we have 40 × 3 + (2 2 + 2 2 ) × 3