key: cord-0883191-273y4eca
authors: Hasaninasab, Mehdi; Khansari, Mohammad
title: Efficient COVID-19 testing via contextual compressive sensing
date: 2021-08-13
journal: Pattern Recognit
DOI: 10.1016/j.patcog.2021.108253
sha: 1f9225ead3edf4c18ff6007eb0914a2990306654
doc_id: 883191
cord_uid: 273y4eca

The COVID-19 pandemic is threatening billions of people's life all over the world. As of March 6, 2021, covid-19 has confirmed in 115,653,459 people worldwide. It has also a devastating effect on businesses and social activities. Since there is still no definite cure for this disease, extensive testing is the most critical issue to determine the trend of illness, appropriate medical treatment, and make social distancing policies. Besides, testing more people in a shorter time helps to contain the contagion. The PCR-based methods are the most popular tests which take about an hour to make the output result. Obviously, it makes the number of tests highly limited and consequently, hurts the efficiency of pandemic control. In this paper, we propose a new approach to identify affected individuals with a considerably reduced No. of tests. Intuitively, saving time and resources is the main advantage of our approach. We use contextual information to make a graph-based model to be used in model-based compressive sensing (CS). Our proposed model makes the testing with fewer tests required compared to traditional testing methods and even group testing. We embed contextual information such as age, underlying disease, symptoms (i.e. cough, fever, fatigue, loss of consciousness), and social contacts into a graph-based model. This model is used in model-based CS to minimize the required test. We take advantage of Discrete Graph Signal Processing on Graph (DSP(G)) to generate the model. Our contextual model makes CS more efficient in both the number of samples and the recovery quality. Moreover, it can be applied in the case that group testing is not applicable due to its severe dependency on sparsity. Experimental results show that the overall testing speed (individuals per test ratio) increases more than 15 times compared to the individual testing with the error of less than 5% which is dramatically lower than that of traditional compressive sensing.

COVID-19 has become the main threat to people's lives all over the world. It also disrupts both business and social activities. Economic depression is another problem that all people are involved in, regardless of affected or nonaffected. In the absence of a vaccine and definite cure, the only remedy is that stop spreading the viruses. Largescale quarantine is not a good solution because of serious, perhaps irreparable consequences in long term including a downturn and social problem. The best policy which minimizes economical damage is to isolate the infected people. So, COVID-19 testing is the primary step to isolate infected people. This testing must be fast enough to cover most people in a short time and consequently isolate the right individuals before time passed. There are mainly two types of test that can detect COVID-19 in individuals: 1) Serological tests, which study the presence of antibodies in the blood of individuals, and 2) The Swap test, which takes material from the cavity between nose and mouth and looks for RNA of an alive virus. Although the serological test has its advantages, the swap test is extremely recommended by the Center for Disease Control and Prevention (CDC) [1] . This class of tests uses Real-Time Polymerase Chain Reaction (RT-PCR) for detecting COVID-19. This method is based on selective DNA strands amplification in each cycle which can make the result in a shorter time and conservation of reagent and testing kits.

Additionally, increasing the range of detection makes it suitable techniques to detect COVID-19. The RT-PCR is done typically in 40 cycles, and in each cycle, three major steps are performed [2] :

High temperature (typically 95°C, the maximum temperature that DNA can withstand) incubation is used to "melt" double-strand DNA to the single stand.

During annealing, the temperature is reduced (usually 5°C below the first step) to let the complementary sequence hybridize.

At 70-72°C, the optimal activity of DNA polymerase is achieved and primer extension occurs.

It is expected to take 8 hours to complete the RT-PCR test for each individual. This is a long time that extremely limits the number of individual testing per day.

Although the RT-PCR test is recognized as the golden standard for determining COVID-19 status in individuals [3] , but it's time-consuming is the critical challenge. It may cause to diagnose from chest CT image as an alternative way [4] .

One of the usual methods to mitigate this limitation is group testing [5] . Group testing is an old technique introduced by Dorfman in 1943. In the group testing, people are grouping and testing is done on the pooled sample.

If the result is negative, it means all group members are disease-free; otherwise, the test must be performed for each group member separately. Although it has advantages compared to sampling techniques (i.e. no data loss), it has two important restrictions: The output of the group testing is a binary value. In other words, it just determines disease affection or non-affection. This is inefficient in the case of COVID-19 where the viral load (i.e. virus concentration) is important to determine the COVD-19 severity among patients [6] and also to study about the antibody response [7] .

The other drawback of group testing is its dependency on the size of groups. Increasing the group size increases the number of required test and fades the efficiency of the group testing [5] . This is extremely critical in case of a large population (e.g. testing .

From the last decade, Compressive sensing (CS) has been a popular technique in signal processing related fields.

The key advantage of CS over traditional sampling is the extremely lower sampling rate required for the signal recovery. CS has two main phases: signal sampling and signal recovery. Contrary to the traditional sampling, the CS samples cannot directly be mapped to any specific data elements, but a linear combination of all data elements. In each CS sample, each individual's data have its own weight which is specified by the sampling matrix. In the recovery phase, the signal elements are recovered from CS samples. From the mathematical point of view, CS signal recovery is going to find a solution (i.e. the original signal) of the underdetermined system of equations that have an infinite number of solutions [8] . It is done by considering the signal's unique features (i.e. sparsity) among other solutions.

In the application of COVID-19 testing, the CS sampling rate is the number of required tests to identify all group members' COVID-19 affection status (i.e. recovered signal). Depending on the recovery strategy chosen, the appropriate sampling rate may be varied [9] . Even though CS efficiency, its implementation has some drawbacks in the real world. In the COVID-19 testing scenario, decreasing the number of required tests while maintaining the required quality is the most critical challenge. Baraniuk introduced the model-based CS to present the unique features concretely in the form of a mathematical model. In model-based CS, in addition to sparsity, the structure of signal values and locations are used as signal features. It has been shown that model-based CS makes it possible to reduce the sampling rate without scarifying robustness and the quality of the output [10] .

To the best of our knowledge, the model-based CS hasn't applied to COVID-19 so far. Moreover, in none of the proposed models, the interdependency between signal elements (i.e. the relation between individual test results in the context of COVID-19) was considered. In this paper, we introduce the Contextual Model-based CS (CMb-CS) that considers the contextual information ( e.g. lifestyle, underlying disease, and social contacts ) to determine the signal elements (i.e. individual test results) interdependency. The model is generated based on the emerging field of Signal Processing on Graph (DSP G ) and exploit its powerful features to show interdependency between signal elements. The aim of our model is to represent the signal characteristic to recover the compressed signal (i.e.

compressed COVID-19 test results) more efficiently. The interdependency between signal element, along with the signal smoothness and sparsity are three factors that are embedded into our model to represent the signal more precisely. Owing the clear image that the model offers from the recovered signal, there are expected that the signal is recovered with much lower number of data. At the final part of this paper, the model efficiency in testing is compared to the existing group testing and the CS based method.

The rest of this paper is organized as follows: Section 2 describes some necessary backgrounds. In Section 3, the proposed model is presented. In section 4 the experimental results are shown and explained, and finally, the conclusion is discussed in section 5.

Group testing was firstly introduced by Dorfman to detect syphilis from blood samples in a large population without the need to test each individual [5] . Instead of testing each individual separately, the population of people is divided into groups of members. In each group, the samples are mixed and the test is done on the pooled samples. The negative result means all group members are free from infection. Otherwise, there is at least one person contaminated. It is proven that in this way, the number of tests per individual is reduced. The efficiency of group testing depends on two key parameters: infection rate ( e.g. probability of infection) and group size ( ). The expected relative cost ( ) given by:

Where is the size of the group (number of members), and is the average infection probability.

Group testing has been applied on COVID-19 testing [11] . The multistage group testing approach was proposed to reduce the number of COVID-19 required tests. In this work, the group testing was done iteratively, and in each iteration, the group size was changed dynamically. It was shown that this approach reduce the above parameter [11] . From (1), it is clear that increasing or may causes , which means that using group testing is ineffective. Therefore, generally group testing is not beneficial in a large population. Moreover, it gives only binary results (infection or non-infection) and the viral load is also unspecified. These shortcomings make the group testing not a suitable method for COVID-19 pandemic conditions (i.e. high population, and fast result required, and high output quality demanding). Another demand is to reduce more number of required tests which consequently, saves quantity of reagents and manpower. Regarding all above together, one of the promising approaches to address these challenges, is compressive sensing which is described in the following section.

CS is an emerging field of research that merges sampling and compression in one step. Contrary to traditional sampling, in the CS, it has been proven that the signal can recover completely at a much lower sampling rate than that in Shannon/Nyquist [8] . Due to the merging of data acquisition and data compression in one single step, there will be no need to access all data elements. One of the key concepts in CS is sparsity which is described in the following.

Imagine the is a vector consist of the N individual test result. Vector x is K-sparse if it has at most K nonzero coefficients (K<<N).

The first step for applying the CS to a signal is the sensing process for which the sensing operator ( ) is used. The compressive measurement is shown as:

(2)

Where the measurement is (sensed) vector. Clearly, each CS sample is a linear combination of all signal elements. The CS technique is a suitable alternative to the Shannon/Nyquist based sampling scheme if the signal is sparse. It can be easily shown that the CS signal recovery problem is solving an underdetermined linear equation system (i.e. the number of unknowns is more than the number of equations). Given that there are an infinite number of solutions, it is necessary to use prior knowledge of the signal to choose the most appropriate solution (x). The sparsity level is one of the key knowledge that is exploited in the traditional CS recovery.

Suppose that is a -sparse signal. The CS recovery process is finding the signal which satisfies two conditions:

first, is -sparse, and second, satisfies (2) . Suppose also that has elements, so finding a -sparse solution simply means searching for a solution among . / different signals ( (3)).

‖ ‖ subject to:

Where ‖ ‖ is the zero norm which is the number of nonzero elements of .

While the CS shows significant improvement from traditional sampling in terms of both sampling rate and output quality, it has some drawbacks. Computing x from . / possible solutions is a time-consuming process, because of the large search space. One possible way to improve the reconstruction process is to limit the search space. It is proven that in this way the reconstruction time is decreased. Furthermore, reconstruction can be done at a much lower sampling rate without sacrificing signal quality [10] (i.e. individual test results).

Limiting search space could be achieved using additional knowledge beyond sparsity. This is the core idea that led to the emergence of the model-based CS, which is a suitable alternative to traditional CS, especially in time and sample constraint applications such as COVID-19 testing. In this way, due to the reduction of the search space too much lower than . /, the recovery process can be done with ever fewer CS measurements. It means that with a lower number of tests, more people COVID-19 status can be identified.

Graph signal processing is an emerging field that tries to express the regular and irregular structured signals [12] .

Social networks and wireless sensor networks are examples of those signals. The signal structure can be represented by a graph form appropriately; the edge weight between two nodes is proportional to the similarity or correlation of them. This unique feature makes the graph a powerful tool to represent the interrelationship between signal elements. In the COVID-19 testing application, the first step is to model the data in the DSP G domain. In the following subsections, the graph signal model and also some of the key concepts in DSP G are discussed.

Consider a graph ( ) where * + is graph node set and * + is a signal value in each node. In graph signal representation, the signal elements can be represented by nodes value and the relations between these elements are represented by graph weight. In the COVID-19 testing scenario, each node corresponds to an individual, is an output vector and its elements is the health status of individual. is the weighted adjacency matrix and is the weight between the nodes i and j (i.e. individual i and j). It is clear that for the undirected graph, is a symmetric matrix (i.e. = ) and =0 if there is no correlation between the nodes i and j. In COVID-19 testing scenario, the weighted adjacency matrix is an appropriate tool to represent the degree of similarity between individuals, which is a part of the proposed model. Figure 1 shows the weighting advantage of DSP G . In the traditional form of signal representation, the relation between signal elements (i.e. red line) is not considered. So all the three signals are the same in traditional signal representation, but they differ in DSP G domain. 

The adjacency matrix of the graph is defined as ( ) in which:

The 1s in the adjacency matrix can be replaced with which typically is the correlation between and .

For graph with adjacency matrix , for the vector the eigenvalue is defined as follows:

The vector is an eigenvector corresponding to .

The graph shift is one of the basic operations in DSP G . Graph shift simply replaces the signal value at the node with the linear combination of signal values in neighborhood nodes.

In traditional signal processing, the Total Variation (TV) is the difference between two consecutive samples [12] . In the case of the COVID-19 testing scenario, TV has a conceptual meaning which shows how the graph signal of test results varies. The TV measures the difference between each node and its neighbors on the graph. In DSP G , TV is generally defined as a similarity between the graph signal and its shifted version. In fact TV is derived from the local variation which somehow represents the counterpart of the signal gradient in traditional signal processing [13] . Graph variation is a highly popular parameter in optimization and learning problems. It characterizes the smoothness of the graph, which is the key parameter to find the discriminative graph signal.

Using the fact that discriminative signals have different values in multiple groups, low smoothness value can be one of the features in that type of graph signal. It can also be used as the key constraint in learning related problems including graph based learning methods [14] . As in [15] , the edge weight smoothness and graph sparsity is considered as two main criteria for dimension reduction and semi-supervise classification. In DSP G , the local variation (gradient) at the n th node is defined as:

Substitute (6) in (7) yields:

Where is a normalized version of A:

Where λ max (A) is the largest eigenvalue of A. Using in (8) guarantees the graph signal is scaled properly compared with the original one.

To make sure that the shifted graph signal is scaled properly compared to the original one, ( ) can be rewritten as:

The magnitude of local variation of the graph in each node is denoted by ‖ ( )‖ [13] . The total variation is naturally the sum of local variations in each node. TV is represented based on ℓp-norm:

For p=2, it is called the quadratic form of graph total variation [16] .

Let us now introduce the CMb-CS for COVID-19 testing. The proposed method can be divided into two main phases: sampling and recovery.

The schematic of the sampling phase is depicted in Figure 2 . In the CS sampling phase, the samples are pooled according to the weights specified in the sampling matrix. The pooled sample is fed in to the RT-PCR machine. So far, it is just similar to the group testing except that the weighting mechanism. The weight values in group testing are the binary values (0 or 1). Although the binary weighting matrix is also applicable in CS-based crowd sensing applications, some of the sensed elements are not contributed in the CS measurement, which may lead to a less reliable result in the recovered signal [17] . . .

Compressed output ( ) .

. . At the recovery phase, the aim is to recover the signal from CS samples. As mentioned before, in the model-based CS, the main goal is to limit the search space by applying a model to find the solution more efficiently. Our As can be seen, the graph is a full mesh graph and each node and the weights correspond to the individuals and the interdependency between individuals' state respectively. As mentioned in the background section, the weights represent the similarity between the contextual information: We now formulate the recovery problem in DSP G . The basis for test results recovery is expressed in (13) . To recover the test results for each individual (x) from the CS measurement data y, x must be a solution for (2) which typically has an infinite number of solutions [9] . To select an individual test results from this infinite number of solutions we must solve a combinatorial multi-parameters optimization problem. In the proposed method, the proper solution (i.e. the individual test result) must minimize the three mentioned limiting factors of the search space and also satisfy (2) at the same time:

Where, ̂ is the recovered signal (i.e. COVID-19 testing result vector), ( ) is the p-norm of graph total variation, ‖ ‖ is the graph ℓ q -norm node value (e.g. node value sparsity), and ‖ ‖ measures the

structural similarity between the recovered signal graph and contextual information graph. This measurement is done simply by making the ℓ r -norm difference between the contextual information graph adjacency matrix ( ) and the recovered data graph adjacency matrix ( ). , and are the weights assigned to each term to adjust their effectiveness on the output result. It is proven that (13) is NP-hard in general and not easy to solve [16] .

Setting proper values for and for each term can express (13) as a convex optimization problem. In the case of graph total variation ( ) although both ℓ 1 norm and ℓ 2 norm (i.e. and respectively) can tend to the convex form, here we use a quadratic norm (i.e. ℓ 2 -norm) instead of an ℓ 1 -norm because; first, it is computationally easier to optimize ℓ 2 -norm than ℓ 1 -norm based graph variation and second, the smooth graph signal can be more efficiently separated from non-smooth parts by ℓ 2 -norm [16] .

However, it is shown that in general ‖ ‖ approaches ‖ ‖ when tends to zero. Therefore, to relax the NP-hard problem, some other values of and are chosen for the second and third terms in (13) , respectively. It is proven that for , the problem again becomes NP-hard. But for =1 it becomes a convex problem [18] . Also, by setting , the problem is still in convex optimization form but it does not provide a reasonable approximation to the original signal [8] . Similarly, r should set to 1. Here, the easier optimization process is the main advantage that motivates us to try to transform the problem into a convex optimization problem since the local minimum is also the global minimum. By setting , , and to 2, 1, and 1 respectively in (13) and (14) is obtained as follows.

As mentioned earlier, in this problem the recovery process is considered in three aspects: graph variation, graph sparsity, and graph structure similarity to contextual information graph. For a simple and efficient solution of (14), a class of decomposition-based methods is considered. These methods are the proper choices especially for solving large-scale statistical tasks. Another advantage is that this class of algorithms is flexible enough to solve multidimension problems [19] .

To solve (14), we apply the Alternating Direction Method of Multipliers (ADMM) in this case, to benefit from the advantages of both decomposability and fast convergence property [19] . Using the ADMM, the minimization problem can be decomposed to the individual sub-problems. As done in a recent work, the problem of minimizing the rank, sparsity, and total variation is solved by ADMM for principal component analysis for compressed sensing of dynamic MRI [20] .

In ADMM, each variable updates iteratively until the convergence takes place (see Appendix 1) . To solve (14) by employing ADMM, we rewrite it as follows.

Where and:

Note that we move from the constraint to the objective function and putting ( ) in the objective function is equivalent to setting to zero. Then we construct the augmented Lagrangian function:

Where is the penalty parameters. We minimize ̂ with respect to each variable individually. For X, we aim to solve:

(See Appendix 1 for details).To minimize (17), we can replace ( ) by (11) and simply take derivation which yields:

Solve for :

For ‖ ‖ , we aim to solve:

Similarity, for minimizing (20) :

Solving for ‖ ‖ : Where ̅̅̅̅̅̅ ‖ ‖ is the axillary variable.

And for :

. ‖ ‖ ‖ ‖ / (23) Also this follows the same as that for ‖ ‖ and :

For simplicity, set :

Solving for :

Where:

Defining and which is cumulative residual at kth iteration:

Equations (17) to (28) are solved iteratively until the two main conditions are met, i.e. the resultant signal satisfies (2) , which is the core equation in compressive sensing, and the , , and remain nearly constant during the iteration. It means:

Where is the threshold for (2). Since ADMM converges slowly in cases that required high accuracy, we chose those thresholds small enough instead of absolute zero to make a tradeoff between signal quality and iteration number.

We obtained the closed-form solution by setting the derivative to zero. We just set to satisfy the constraints.

For the final implementation, see Algorithm1.

Input : , , ,

While the stopping criterion is not satisfied 

The proposed model is applied on data with the clinical characteristic of COVID-19 infected individuals in Iran [21] .

This dataset contains 1024 individuals' data which covers a wide variety of ages (15-85 years) . In this dataset, the affection rate is 0.46 (i.e. 472 individuals were affected). In the case of group testing, the expected relative cost ( ) in (1) is greater than one and so, the group testing is totally useless. The traditional CS and the proposed CMb-CS is The results clearly show that applying our proposed model, with just 150 tests, can identify all 1024 individuals' affection status correctly. While in traditional CS a similar result still can't be achieved even by increasing fourfold the test numbers. It is due to the fact that traditional CS mainly relies on sparsity, which is not a highly descriminative feature in this case. In addition to the sparsity, our proposed model considers some real features (i.e. graph smoothness and graph similarity to the reference graph model). In Table 1, along with MSE, false positive  and false negative results are represented at the various test numbers ranging from 50 to 200. As the table shows, our proposed method outperforms traditional CS in all number of tests.

The effect of each contextual information item is represented in Figure 6 . It is clear that "age" and "loss of consciousness" have a major and minor effect on the result respectively. Comparing to the Traditional CS, it is clear that even if the CMb-CS uses only the least effective contextual information (i.e. "loss of consciousness"), it still outperforms traditional CS. Figure 6 also reveals that the "age" and "fever" have considerable effect compare to the other contextual information and therefore, are the most dominant factors for identifying COVID-19.

Regarding other contextual information, they are more separable from the sample rate of 260 onward. sparsity, graph smoothness, and graph similarity respectively. As can be seen from Figure 7 , changing the sparsity component weight results in lower change in the output result, which implies that in this scenario, sparsity is not a critical parameter. On the other hand, it seems that the changing graph similarity component weight has a considerable change in the output result (See Figure 9 ), which means that this is the most effective component in this experiment.

The effect of eliminating each model component is also represented in Figure 10 . conclusion In this paper, the contextual model-based CS is introduced to decrease the number of required tests for detection . Despite the popularity of group testing, it was shown that this approach is not applicable in most realworld scenarios because of its strong dependency on data sparsity. Although the CS outperforms group testing, it also suffers from sparsity dependency which limits its efficiency in real scenarios. 

Consider the following convex optimization problem:

Where , and is convex. The Lagrangian for solving this problem is:

Generally, the optimization problem can be examined from two aspects: primal problem and dual problem. The solution for the dual problem provides a lower bound for the primal solution. The main advantage of this approach is that we can solve the problem regardless of how complex the primal problem is. The dual problem is easy to solve due to its concave minimization nature. Using dual function would also lead to an efficient and reliable solution.

the dual problem is defined as :

Where is a dual function of [19] . Thus, the dual problem is:

We can find a solution for the primal problem ( ) from dual problem solution ( ) as follows:

( )

Several methods can be used to solve the above problem. The first one is the dual ascent that solves it using gradient ascent. The dual ascent method is implemented as follows:

Where the step size and k is the iteration counter. and are the primal and dual variable, respectively.

One of the main advantages of dual ascent is that it can be implemented in a decentralized form. In other words, is separable in the form below:

Where ( ) and are subvector of x. The matrix also can be decomposed as follows:

Hence, ∑ . The Lagrangian problem can be written as:

In this way, the a-minimization step can be divided into N separate problems and solved in parallel as follows Despite the merit of the dual ascent method, it has some limitations. The main disadvantage is that if f is a nonzero affine function of x, then L is unbounded below in x and so the x-update fails.

The Augmented Lagrangian methods were introduced to add robustness to the dual ascent method and also reduce the limitation of f. The Augmented Lagrangian is:

Where is the penalty parameter. Applying The Augmented Lagrangian method yields the following algorithm:

This form is known as the method of multipliers and it has been proven that it covers more general conditions. But this advantage comes at a cost. The is not separable even if the is. Hence, the minimization step cannot be carried out in parallel to x.

ADMM is the algorithm that benefits from both decomposability of dual ascent and good convergence of the method of multipliers. ADMM solves the problem of the form:

Where and are both convex. As can be seen, in this form the variable has been split into two parts ( and ) and also the objective function is separable. The augmented Lagrangian is

The ADMM algorithm finds the optimal solution via iterations:

In ADMM, x and z are updated sequentially. The minimization of each variable in separate steps is what allows decomposition when f or g is separable.

ADMM scaled form ADMM can be rewritten differently by combining the two latter terms in (38) and scaling the dual variable( ):

Where is a scaled dual variable and . We can express the ascaled form of ADMM as follows:

. Where is the cumulative residual:

☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Overview of COVID-19 testing and implications for otolaryngologists

Real-time PCR handbook

Synergistic learning of lung lobe segmentation and hierarchical multi-instance classification for automated severity assessment of COVID-19 in CT images

Momentum contrastive learning for few-shot COVID-19 diagnosis from chest CT images

The detection of defective members of large populations

Nasopharyngeal viral load predicts hypoxemia and disease outcome in admitted COVID-19 patients

Kinetics of viral load and antibody response in relation to COVID-19 severity

An introduction to compressive sensing

A mathematical introduction to compressive sensing

Model-based compressive sensing

Multi-stage group testing improves efficiency of large-scale COVID-19 screening

Discrete signal processing on graphs: Frequency analysis

Local Smoothness of Graph Signals

A feature extraction model based on discriminative graph signals

Sparse graphs with smoothness constraints: Application to dimensionality reduction and semi-supervised classification

Signal recovery on graphs: Variation minimization

Binary sparse signal recovery algorithms based on logic observation

A note on the complexity of L p minimization

Distributed optimization and statistical learning via the alternating direction method of multipliers

Smooth robust tensor principal component analysis for compressed sensing of dynamic MRI

Data on epidemiological and clinical characteristics of coronavirus-infected (COVID-19) individuals in Iran

Age-dependent effects in the transmission and control of COVID-19 epidemics

Global Health, a data science initiative