key: cord-0666702-djbp78lf authors: Petersen, Hendrik Bernd; Bah, Bubacarr; Jung, Peter title: Practical High-Throughput, Non-Adaptive and Noise-Robust SARS-CoV-2 Testing date: 2020-07-17 journal: nan DOI: nan sha: f61f75a7365139fa5874feee53983b8796f035af doc_id: 666702 cord_uid: djbp78lf We propose a compressed sensing-based testing approach with a practical measurement design and a tuning-free and noise-robust algorithm for detecting infected persons. Compressed sensing results can be used to provably detect a small number of infected persons among a possibly large number of people. There are several advantages of this method compared to classical group testing. Firstly, it is non-adaptive and thus possibly faster to perform than adaptive methods which is crucial in exponentially growing pandemic phases. Secondly, due to nonnegativity of measurements and an appropriate noise model, the compressed sensing problem can be solved with the non-negative least absolute deviation regression (NNLAD) algorithm. This convex tuning-free program requires the same number of tests as current state of the art group testing methods. Empirically it performs significantly better than theoretically guaranteed, and thus the high-throughput, reducing the number of tests to a fraction compared to other methods. Further, numerical evidence suggests that our method can correct sparsely occurring errors. During the early stages of a pandemic an exponential growth in infected individuals [1] might result in a shortage of tests for an infection. In later stages of a pandemic regular testing of large groups of mostly healthy people is a common policy to prevent further outbreaks [2] [3] . This results in a strong desire to reduce the number of tests required to identify infected elements among a group of individuals. This approach is often called pooling and has been intensively investigated for COVID-19, see for example [4, 5, 6] . The mathematical field group testing is concerned with this problem and has therefore also gained interest recently again. It was invented by Dorfman [7] , however for a general overview we refer to [8] . The methods from classical group testing often suffer from several drawbacks including slowness due to adaptivity of the tests [9] and sensitivity to errors [10] . Compressed sensing is a mathematical field concerned with recovering a vector with many zero entries from as few as possible non-adaptive linear measurements. For a general overview we refer to [11] . Compressed sensing has achieved several theoretic goals including achieving the optimal number of measurements [11, Chapter 13] , independence of the measurements from each other, robustness to noise [11, Chapter 4] and sometimes even error correcting properties [12] . We propose to use compressed sensing instead of classical group testing for viral detection since it tackles some weak points in group testing mentioned above. Our contribution in this survey paper is as follows. We explain how the viral detection problem can be modeled as a compressed sensing problem. We present one recent result from compressed sensing to derive a complete testing procedure for viral detection. Our method will require the same number of tests as state of the art group testing methods in several examples. However, our approach is non-adaptive and has rigorous recovery bounds showing robustness to noise. Numerically we see that our method has empirically further advantages, including error correcting properties. Compared to other papers concerned about using compressed sensing for group testing [13] [14], we use a different compressed sensing decoder that is more robust against outliers [15] . The rest of the paper is organised as follows. In Section 2 we model the viral detection problem as a compressed sensing problem. Precisely, we present in Subsection 2.1 the pooling design, in Subsection 2.2 the underlying matrices for the pooling design, and in Subsection 2.3 the identification of infected individuals. We discuss the advantages of our approach and give a numerical advantage of our approach in Section 3. For completeness we present the key mathematical details supporting the approach we proposed in Section 4. Suppose we want to find S individuals infected with a virus among N individuals. According to the information theoretic lower bound we require at least ln N S tests to find the infected individuals. The binary splitting algorithm finds the S infected individuals with a number of tests in the order of S log 2 (N/S). However, since the tests are adaptive, each subsequent test is designed depending on the outcome of previous tests, and thus each test has to wait for the result of previous once. If the time to perform a test is large, it might be desirable to perform multiple tests at once. In such cases non-adaptive methods are preferable. There are many deterministic non-adaptive methods for group testing. Most of these prove their results using disjunct matrices. For instance, in [16] a non-adaptive group testing method is presented whose number of tests is in the order of S 2 log 2 (N )/log 2 (S). We propose to use compressed sensing with additional nonnegativity constraints to solve this problem. We collect specimens from N individuals arranged into a vector x ∈ R N , and the amount of viruses in the specimen of the n-th individual is denoted by the non-negative quantity x n ∈ N ∪ {0} ⊂ R. We assume that the viruses are evenly distributed in each specimen, meaning that if we take α ∈ [0, 1] of the volume of the specimen of the n-th individual, it will contain roughly αx n viruses. We denote the number of non-zeros in a vector z by z 0 . Since, S individuals are infected we have x 0 = S. Let the sample of the m-th test contains a fraction of A m,n of the amount of specimen of the n-th individual. The sample of the m-th test thus contains, up to rounding errors, the amount of viruses n∈[N ] A m,n x n = (Ax) m , where A is an M × N matrix with entries A m,n ∈ [0, 1] with column sums of at most one. The amount of tests is thus M . Often one makes the assumption that x is a random vector with iid. elements, for instance a Poisson random variable multiplied by a Bernoulli random variable with parameter p = S N . Instead we take it here as deterministic unknown. A quantitative real-time reverse transcription polymerase chain reaction (qPCR) [17] can be used to generate an estimate y m of the amount of viruses in the m-th test (Ax) m . This procedure is not accurate and errors e m := y m − (Ax) m might occur. We try to recover the amount of viruses in the specimen of the individuals according to with a possibly small x 0 and a as small as possible M . This is exactly a compressed sensing problem where, due to the nature of the problem, x is non-negative and exactly S-sparse (compressed sensing guarantees often also work for compressible vectors which are well approximated by sparse vectors). The theory of compressed sensing states that there exists indeed matrices and efficient decoders that allow recovery of x if M is in the order of S ln e N S [11] . It remains to choose a suitable matrix and a decoder. We begin with the sensing matrix. The theory of compressed sensing states that there exists indeed matrices that allow recovery of x from (1) if M is in the order of S ln (N/S), referred to as the optimal scaling. This testing design described above is a special compressed sensing task where the elements of the matrix A and the unknown sparse vector x are nonnegative. For non-negative A standard proof techniques like the restricted isometry property does not succeed out of the box in this regime and instead one has to use tools like the null space property defined in Section 4. For example, it is known that matrices with independent and uniformly on {0, 1} distributed entries achieve the optimal scaling [18] . However, deterministic, i.e. non-random, construction of matrices with a null space property is a delicate, open problem. To the best of the authors knowledge the best deterministic construction of binary compressed sensing matrices is [19, Theorem 3.5] , which generates expander graphs that are sufficient for compressed sensing according to [11, Chapter 13] . In particular for any α > 0 there exists a constant C α such that the matrix has only M ≤ C α S 1+α (log 2 (N ) log 2 (S)) 2+ 2 α rows. The downside of this result is that the constant C α is rather large for small α and thus other explicit constructions often achieve a smaller number of tests, since the dimensions for viral detection might be rather limited. We therefore, propose the binary matrix and the S-disjunct matrix constructions in [20] and [21] respectively, in which the number of measurements M is in the order of S √ N . Actually, we established the equivalence of the two constructions in Section 4. Such S-disjunct matrices do not achieve the optimal rate for fixed S as N → ∞ according to [22] . However, for viral detection we have a fixed prevalence p = S N in mind. The prevalence might be larger in the beginning of a pandemic and smaller when healthcare professionals are testing regularly and asymptotically, but for each of these applications we consider a fixed p = S N . In this scenario our result achieves the rate M = pN 3 2 + N 1 2 , which can outperform results that achieve optimality according to [22] . For instance the construction of [16, Corollary 1] achieves a number of measurements that is less than . For a suitably chosen N , this upper bound and other constructions can be outperformed by the construction we propose. A different idea is to generate matrices at random and test them for a null space property until we find a good matrix. Since many random matrices obey the null space property in the optimal order [11, Chapter 13] , this approach should succeed. However such approaches fail in larger dimensions, since testing for a null space property or some related concepts is in general not possible in polynomial time [23] . The classical decoding procedure for disjunct matrices iterates over all n ∈ [N ] and, if for a fixed n and for all m with A m,n = 0 the m-th test is positive, it declares the n-th individual as infected. This process is simple and, if A is S-disjunct, there is no noise and there are no more than S infected individuals, it is guaranteed to find exactly all infected individuals [21] . However, every false negative test will result in at least one individual that is falsely flagged as not infected. Thus, this decoding procedure is incredibly sensitive to noise. There are extensions to these matrices which can tolerate a fixed number of errors under restrictive conditions [16] . Compressed sensing on the other hand can not only detect infected individuals but also estimate the viral load, which may have further benefits to the medical practitioners. The noise vector in (1) is non-zero in general, since the estimate is affected by some errors including, rounding errors and inaccuracy of the qPCR. In [14] the noise vector is modeled as a heavy tailed random variable depending on the unknown quantity Ax. Therefore, it is difficult to apply recovery methods from compressed sensing for independent additive noise out of the box. Parameter tuning, using e.g. cross validation, is often crucial at this point. Interestingly, non-negativity helps with such noise models. Combining the heavy-tailed noise model with non-negativity of data, we arrived at a conclusion that the best is to use the parameter tuning free Non-negative Absolute Deviation Regression (NNLAD) for recovery which is any minimizer In [15] we showed that for certain matrices A the convex NNLAD approach indeed is sparsity promoting and allows for an estimate of the form x − x # 1 ≤ C e 1 , for some constant C independent on e, x, y and x # . Finally, given a certain threshold ǫ one declares the individual n to be infected if x # n > ǫ. If the noise level e 1 is small enough, this method will guarantee that these are exactly the infected individuals. Thus, in this sense compressed sensing gives a non-adaptive testing procedure which provably finds S infected individuals among N with the scaling-optimal number of tests, and small errors would have no effect on the test result. The guarantees in [15] follow from a self-regularization feature in the non-negative case which has been worked out already for other cases, like the non-negative least squares in [18] . We compare our method to other methods proposed for group testing. Consider the non-adaptive method from [16] . Their number of tests is less than p 2 N 2 log 2 (N ) log 2 (N )+log 2 (p) tests. For a prevalence of 0.01 at 900 individuals their method uses at most 0.5573 tests per individual while ours uses 0.3333 tests per individual at a pool size of 30. For a prevalence of 0.001 at 10000 individuals their method uses at most 0.08 tests per individual, while ours uses 0.11 tests per individual at a pool size of 30. In many scenarios our method requires a similar amount of tests as other methods from non-adaptive group testing and sometimes even requires less tests. Adaptive methods outperform non-adaptive methods in general. Consider the adaptive method from [7] . For a prevalence of 0.01 at 900 individuals the approach of [7] expects 0.2 tests per individual at a pool size of 11 while our method uses 0.3333 tests per individual at a pool size of 30. For a prevalence of 0.001 at 10000 individuals the approach of [7] expects 0.0628 tests per individual at a pool size of 32 while our method uses 0.11 tests per individual at a pool size of 100. The method of [6] is even better than [7] . At a prevalence of 0.01 their method only requires 0.12 tests per individual at a pool size of 35, and at a prevalence of 0.001 only requires 0.018 tests per person at a pool size of 350. Consider the adaptive method from [24] . For a prevalence of 0.01 and at 5000 individuals their method uses 0.3 tests per person. Our method can handle a prevalence of 0.01 at 900 individuals with the same number of tests per person. Many adaptive group testing methods require half as much tests as our method, but ours the advantage that it only requires one stage and can thus be performed faster. We shortly demonstrate an advantage of the NNLAD approach to other compressed sensing decoders in a short experiment with synthetic data. Following Theorem 4.5 for Q = 31 and S = 7, we are guaranteed to detect up to 7 infected people among 961 by using 248 tests. However, empirically the identification will succeed even if more than 7 individuals are infected and even if multiple measurements are corrupted. In Figure (1) we vary the prevalence p = and plot the probability that the NNLAD estimator x # is sufficiently close to the true signal x in the ℓ 1norm. We see that as guaranteed by Theorem 4.5 the recovery succeeds for p = x 0 N ≤ 7 31 2 ≈ 0.0073 and p e = 0. But empirically the recovery also succeeds for p ≤ 0.08 and p e = 0, i.e. ten times higher than guaranteed. This suggests that S and thus M could be reduced and still recover whenever p ≤ 0.02 for instance, which might be sufficient for a prevalence of p = 0.01. Further the NNLAD seems to successfully recover even in the presence of sparse noise, i.e. when p e is small. This is not suprising as the NNLAD minimizes the ℓ 1 norm of all possible noises and ℓ 1 -minimization is known to be sparsity promoting [11, Chapter 4] . Recovery seems to succeeds whenever p e ≤ 0.06 and 4 3 p e + p ≤ 0.08. This error correcting property gives the NNLAD decoding approach advantages over other compressed sensing decoders in the presence of heavy outliers. However, the error correcting properties cannot be guaranteed uniformly for all e with p e ≤ 0.06. If the noise components are active exactly on the support of a column of A where the signal x is non-vanishing, the measurement might corresponds to a different signal with the same support. Hence, recovery might fail as soon as e 0 is at least as large as the number of non-zero entries in a column of A, in this case S + 1 = 8. However, there are M S+1 = 248 8 ≈ 3.16 · 10 14 different possible supports of a S + 1-sparse noise, but only N = 961 of those appear as columns of the matrix. Thus, if the noise support is drawn uniformly at random, such an event is highly unlikely. Thus, we see that the NNLAD is empirically correcting more errors than expected. We have explained how compressed sensing can be used to solve the viral detection problem. It generates a non-adaptive testing procedure. We have seen that certain design matrices from non-adaptive group testing are already suitable for compressed sensing based viral detection. Further, we have presented a construction of design matrices that can be used for classical non-adaptive group testing and compressed sensing based viral detection. The construction requires roughly as much tests as other methods from non-adaptive group testing and can possibly outperform those. Adaptive group testing methods still require less tests but can only be performed sequentially which is a critical problem in a pandemic with an exponential growth rate. We have proposed to use the NNLAD as a compressed sensing decoder, since, compared to other compressed sensing decoders, it is robust against heavy tailed noise and does not requires knowledge of noise level. In particular, it also obeys certain error correcting properties. Lastly our method also computes the viral load of infected individuals. Classical, deterministic, non-adaptive group testing often makes use of so called disjunct matrices. holds true. Then, A is called S-disjunct. The classical decoding procedure for a disjunct matrices iterates over all n ∈ [N ] and, if for all m ∈ A n the m-th test is positive, it declares the n-th individual as infected. This process is simple and, if A is S-disjunct and there are no more than S infected individuals, it is guaranteed to find exactly all infected individuals [21] . However, every false negative test will result in at least one individual that is falsely flagged as not infected. Thus, this decoding procedure is incredibly sensitive to noises. There are extensions to these matrices which can tolerate a fixed number of errors by forcing the set A n \ n ′ ∈T A n ′ not only to be non-empty but also sufficiently large [16] . Compressed sensing on the other hand makes use of matrices that have a null space property. holds true. Then we say A has the robust null space property of order S with constants ρ and τ . If A has the robust null space property of order S any vector x with at most S non-zero components can be recovered from Ax by solving a linear program and, in the presence of noise, the estimation error is bounded by a constant times the norm of the noise [11, Chapter 4] . It is quite suprising that there is a common way to generate S-disjunct matrices and matrices which have the robust null space property of order S. This was proven in [21, Equation (4) ] and [20] . Proof. Note that λ is also the maximal number of common ones of two distinct columns of A. According to [21, Equation (4) where v ′ is the projection of v onto the null space of A and v ′′ is orthogonal to the null space. Since the Moore-Penrose inverse composed with the matrix itself is the identity on the orthogonal complement of the null space, we get By [25, Theorem 2] or [20, Theorem 7] for any n ∈ [N ] we have |v ′ n | ≤ λ 2D v ′ 1 . It follows that where the last inequality follows from (4). We set α := 2D λ and β := λ 2D + 1 A † 1→1 . Using (6) for all n ∈ [N ] and v ∈ R N and noting that S = D−1 λ < α 2 holds true, allows us to apply [20, Lemma 2] . This yields that A has the robust null space property of order S with constants S α−S = S 2D λ −S = ρ, and This completes the proof. This yields that a lot of test schemes for non-adaptive group testing can already be used with the compressed sensing decoding strategy presented in Section 2. Using this method we can not improve the amount of tests required directly, but, the compressed sensing decoder we propose is significantly more robust to noise. We will see empirically that the compressed sensing decoder will correctly identify all infected individuals even in the presence of one false negative test, unlike the classical group testing decoder for disjunct matrices. Note that under the assumptions of the theorem the quantity λ D is one over the mutual coherence of the matrix (after suitable normalization). The mutual coherence is a general tool in compressed sensing and it is known that the order of the null space property has to be at least in the order of one over the mutual coherence [26, Equation (9) ]. In the special case of binary matrices this result can be refined to account for a better constant using [20] . We introduce the result [15, Theorem 3.4] with some fixed parameters that we will combine with some binary matrices from [20, Page 3015] that have a high number of ones per column but a small inner product between distinct columns. (1) A has the robust null space property of order S with constants ρ ∈ [0, 1) and τ ∈ [0, ∞). (2) There exists some vector t ∈ R M such that A T t > 0 and κρ < 1, where κ := Then, the for all x ∈ R N , x ≥ 0, x 0 ≤ S, e ∈ R M and y := Ax + e we have that any minimizer x # of We now have the necessary tools to proof the main result. For all x ∈ R N + with x 0 ≤ S, y ∈ R M and e := y − Ax any minimizer x # of min z∈R N + Az − y 1 obeys where A † is the Moore-Penrose inverse of A. Proof of Theorem 4.5. Statement 1 is clear. Since P is a permutation matrix, it has exactly one non-zero entry per row and column. Since A has exactly (S + 1) of those blocks per column and Q of those blocks per row, statement 3 follows. Statement 2 follows from 1 and 3. Consider the matrix B := (S + 1) A with columns B n whose entries are either zero or one by construction. It has exactly D = S+1 ones per column. By [20, Page 3015] two distinct columns of B have a scalar product of at most one. Hence λ = max n,n ′ ∈[Q 2 ]:n =n ′ B n , B n ′ = 1, which is also the maximal number of common ones in two different columns. By Theorem 4.3 B is S-disjunct and has the robust null space property of order S with constants ρ ′ = S which yields the claim. Suppose we fix a threshold ǫ > 0 that identifies infected persons, meaning that the n-th person is infected if and only if more than ǫ viruses are contained in the specimen of the n-th person, i.e. if and only if x n > ǫ. In this case we can identify the infected persons even in the presence of small noise. If e 1 < ǫ 4 S + 1 + S (2S + 3) A † 1→1 −1 , we get Since x n is either greater than ǫ or zero, we can deduce that x n > ǫ happens if and only if x # n > ǫ 2 . After the tests we could declare that a person is infected if x # n > ǫ 2 , and healthy if this is not fulfilled. This method would still detect the infected individuals in the presence of small noise. Preliminary estimation of the basic reproduction number of novel coronavirus (2019-nCoV) in China, from 2019 to 2020: A data-driven analysis in the early phase of the outbreak Testing Guidelines for Nursing Homes Prvention und Management von COVID-19 in Alten-und Pflegeeinrichtungen und Einrichtungen fr Menschen mit Beeintrchtigungen und Behinderungen FACT-Frankfurt adjusted COVID-19 testing-a novel method enables high-throughput SARS-CoV-2 screening without loss of sensitivity Evaluation of COVID-19 RT-qPCR test in multi-sample pools A strategy for finding people infected with SARS-CoV-2: optimizing pooled testing at low prevalence The Detection of Defective Members of Large Populations Group Testing: An Information Theory Perspective A Comparative Survey of Non-Adaptive Pooling Designs Pooling designs and nonadaptive group testing: important tools for DNA sequencing A Mathematical Introduction to Compressive Sensing. Birkhäuser Basel Exact signal recovery from sparsely corrupted measurements through the Pursuit of Justice Low-Cost and High-Throughput Testing of COVID-19 Viruses and Antibodies via Compressed Sensing: System Concepts and Computational Experiments A Compressed Sensing Approach to Group-testing for COVID-19 Detection Efficient Noise-Blind ℓ 1 -Regression of Nonnegative Compressible Signals New Construction of Error-Tolerant Pooling Designs Quantitative real-time RT-PCR data analysis: current concepts and the novel gene expressions C T difference formula Robust Nonnegative Sparse Recovery and the Nullspace Property Unbalanced Expanders and Randomness Extractors from Parvaresh-Vardy Codes Compressed Sensing Using Binary Matrices of Nearly Optimal Dimensions Nonrandom binary superimposed codes Bounds on the Length of Disjunctive Codes The Computational Complexity of the Restricted Isometry Property, the Nullspace Property, and Related Concepts in Compressed Sensing Noisy Pooled PCR for Virus Testing Reconstruction guarantee analysis of binary measurement matrices based on girth Sparse representations in unions of bases The work was supported by DAAD grant 57417688. BB has been supported by BMBF through the German Research Chair at AIMS, administered by the Humboldt Foundation.