Message Boards

GROUPS:

Staff Picks Biological Sciences Data Science Medical Sciences Graphics and Visualization Wolfram Language Modeling

Diego Zviovich

Epidemiological Model for repetitive rapid testing for COVID-19

Diego Zviovich

Posted 11 months ago

4976 Views

2 Replies

9 Total Likes

Follow this post

MODERATOR NOTE: coronavirus resources & updates: https://wolfr.am/coronavirus

Epidemiological Model for repetitive rapid testing for COVID19

Diego Zviovich
September 13, 2020

Introduction

Nassim Nicholas Taleb posted a very interesting tweet a couple of days ago.

Given that rapid inexpensive tests, like

Abbot'sBinaxNOW^TM

, which is estimated to cost 5 USD/EA, will be available soon, it is important to explore the ramification of overactive testing as an alternative to stringent lock-downs, as observed in countries like Argentina, or in several states across the US.

Types of tests available

Real-Time Reverse Transcriptase (RT)–PCR Test

The RT-PCR test detects the SARS-CoV-2 virus in upper and lower respiratory specimens. Tests results take approximately 3-5 days to be reported back to the patients.

Searching for information on the type I error (false positive), the median error reported between 2.3% and 5%.

Regarding type II error (false negative), the literature reports a 67% median 4 days after the infection, 38% false positive at the beginning of the symptoms with the low 20% reached 3 days after the onset (~ 8 days after the onset of symptoms).

I was unable to find any information of either false negatives on asymptomatic patients.

In[]:=

solFit=FindFit[{{4,0.67},{5,0.38},{8,0.2},{10,0.2},{20,.2}},dExp[-a(t-b)]+c,{a,b,c,d},t]

Out[]=

{a0.956174,b1.68847,c0.196381,d4.32205}

In[]:=

f[t_]:=Evaluate[dExp[-a(t-b)]+c/.solFit]

In[]:=

Plot[dExp[-a(t-b)]+c/.solFit,{t,2,10},PlotRange{0,1}]

Out[]=

In[]:=

Clear@g

In[]:=

g[t_]:=dExp[-a(t-b)]+c/.solFit

In[]:=

cut=t/.First@NSolve[g[t]1.,t,Reals]

Out[]=

3.44794

We’ll calculate the false negative probability for the slow test with the function
fs.

In[]:=

fns[t_]:=Piecewise[{{1.,t≤cut},{g[t],t>cut}}]

In[]:=

Plot[fns[t],{t,0,10},PlotRange{0,1},PlotLabel"RT-PCR Test -False Negative Ratio",PlotTheme"Detailed",FrameLabel{"Time","False Negative Ratio"}]

Out[]=

fns(t)

Rapid test via Antigen qualitative detection device

On August 26th, the US FDA issues and Emergency Use Authorization for a $5, 15 minute, Covid-19 antigen test, that requires no instrumentation, with demonstrated sensitivity of 97.1% and specificity of 98.5% in clinical study. Paired with a no-charge complementary phone app, which allows people to display their BinaxNOW test results when asked by organizations where people gather, such as workplaces and schools.

In the literature, sensitivity is also referred as positive agreement while specificity is known as negative agreement.

This will allow us to calculate the false negative error in the rapid test.

In[]:=

Out[]=

The false negative for the quick test will be defined by the function fq.

In[]:=

fnq[t_,sensitivity_:0.971]:=(1-fns[t])(1-sensitivity)+fns[t]

In[]:=

Plot[{fns[t],fnq[t]},{t,0,10},PlotRange{0,1},PlotLabel"False vs Slow -False Negative Ratio",PlotTheme"Detailed",FrameLabel{"Time","False Negative Ratio"}]

Out[]=

	fns(t)
	fnq(t)

False positive for the RT-PCR test is estimated to be 7.3%

. Given that the specificity is 98.5%

In[]:=

fps=0.073;fpq=(1-fps)(1-0.985)+fps

Out[]=

0.086905

Modeling the different scenarios with Agent Based Modelling

Model Description

The availability of a rapid test that can be performed by any individual without the intervention of a health worker opens the possibility to allow for events or locations where people gather. By performing the test, and certifying through a mobile application that the user has no detectable antigens, the person might then be allow to engage with others in a close space event.

The scenario described above would require the following protocol:

◼

An individual would submit to a rapid test to determine the presence of antigens on a periodic basis.

◼

In the case that the test returns a positive result, the individual would move to isolation and will perform a detailed RT-PCR slow test to determine if the individual remains infectious or not (the rapid tests shows only if the person was exposed to the virus).

◼

People fully recovered will not need to be submitted for further testing

◼

It is expected that a certain ratio of individuals will break the protocol outlined above: we can expect individuals tat will avoid being tested. Also, a proportion of individuals will break the quarantine.

Compartmental Model

Out[]=

Contents cannot be rendered at this time; please try again later

In[]:=

compartments={{"Comparment","Description"},{"S","Susceptibles"},{"SQ","Susceptibles Quarantined"},{"E","Exposed"},{"EQ","Exposed Quarantined"},{"IA","Infected Asymptomatic"},{"IAQ","Infected Asymptomatic Quanrantined"},{"IS","Infected Symptomatic"},{"ISQ","Infected Symptomatic Quarantined"},{"R","Recovered"},{"TS","RT-PCR Tests"},{"TF","Rapid Tests"}};

In[]:=

rCompartments=Rule@@@Rest[compartments];

In[]:=

ResourceFunction["GridTableForm"][Rest@compartments,TableHeadingsFirst[compartments]]

Out[]=

#	Comparment	Description
1	S	Susceptibles
2	SQ	Susceptibles Quarantined
3	E	Exposed
4	EQ	Exposed Quarantined
5	IA	Infected Asymptomatic
6	IAQ	Infected Asymptomatic Quanrantined
7	IS	Infected Symptomatic
8	ISQ	Infected Symptomatic Quarantined
9	R	Recovered
10	TS	RT-PCR Tests
11	TF	Rapid Tests

In[]:=

parameters={{"Parameter","Description"},{"β","Contact Rate"},{"r","Quarantine Attenuation Factor"},{"lp","Latency Period"},{"ip","Infectious Period"},{"ar","Asymptomatic Ratio"},{"trq","Rapid Testing Execution Ratio"},{"trsS","Slow Testing Execution Ratio symptomatic (use when trq=0)."},{"trsA","Slow Testing Execution Ratio asymptomatic (use when trq=0)."},{"qr","Quarantine Compliance Ratio"},{"fps","False Positive Ratio RT-PCR Test"},{"fpq","False Positive Ratio Rapit Test"},{"fns","False Negative Ratio RT-PCR Test"},{"fnq","False Negative Ratio Rapit Test"},{"cs","Unit Cost per RT-PCR Test"},{"cf","Unit Cost per Rapid Test"},{"g","Graph representing the contact network"},{"tMax","Maximum time frame"},{"time2SlowTest","Time between quick and slow test"},{"timeSlowTestWait","Wait time to obtain slow test results"}};

In[]:=

ResourceFunction["GridTableForm"][Rest@parameters,TableHeadingsFirst[parameters]]

Out[]=

#	Parameter	Description
1	β	Contact Rate
2	r	Quarantine Attenuation Factor
3	lp	Latency Period
4	ip	Infectious Period
5	ar	Asymptomatic Ratio
6	trq	Rapid Testing Execution Ratio
7	trsS	Slow Testing Execution Ratio symptomatic (use when trq=0).
8	trsA	Slow Testing Execution Ratio asymptomatic (use when trq=0).
9	qr	Quarantine Compliance Ratio
10	fps	False Positive Ratio RT-PCR Test
11	fpq	False Positive Ratio Rapit Test
12	fns	False Negative Ratio RT-PCR Test
13	fnq	False Negative Ratio Rapit Test
14	cs	Unit Cost per RT-PCR Test
15	cf	Unit Cost per Rapid Test
16	g	Graph representing the contact network
17	tMax	Maximum time frame
18	time2SlowTest	Time between quick and slow test
19	timeSlowTestWait	Wait time to obtain slow test results

Parameters

In[]:=

aStates=AssociationThread[compartments[[2;;-3,1]],Range[0,Length@compartments-4]]

Out[]=

S0,SQ1,E2,EQ3,IA4,IAQ5,IS6,ISQ7,R8

In[]:=

aColorRules=Association[0Blue,1Darker@Blue,2Yellow,3Darker@Yellow,4Brown,5Darker@Brown,6Red,7Darker@Red,8Purple];

In[]:=

aStatesReverse=AssociationMap[Reverse,aStates]

Out[]=

0S,1SQ,2E,3EQ,4IA,5IAQ,6IS,7ISQ,8R

In[]:=

lAggCompartments={"Susceptibles","Exposed","Symptomatic","Asymptomatic","Infected","Recovered"};

In[]:=

lTests={"SlowTest","FastTest"};

In[]:=

aParameters=AssociationThread[parameters[[2;;,1]],{2.2/100,.2,LogNormalDistribution[Log[3],Log[1.35]],LogNormalDistribution[Log[10],Log[1.5]],4./5.,.0,0.8,.2,.75,fps,fpq,fns,fnq,5.,100.,RandomGraph[WattsStrogatzGraphDistribution[1000,0.04,9]],360,4,3}];

Functions

Algorithm Description

◼

Identify infected cells

◼

Identify susceptible neighbors

◼

Determine if contagion succeeds.

◼

Perform rapid tests > Change agent status

◼

Perform slow tests > Change agent status.

◼

Determine if still infected.

Experiments

Experiment 1 - No rapid test availability

No rapid tests available, 20% asymptomatic get tested with RT-PCR, 80% of the symptomatic, 75% respect quarantine.

In[]:=

ResourceFunction["GridTableForm"][List@@@Normal[aParameters],TableHeadings{"Parameter","Value"}]

Out[]=

#	Parameter	Value
1	β	0.022
2	r	0.2
3	lp	LogNormalDistribution[Log[3],0.300105]
4	ip	LogNormalDistribution[Log[10],0.405465]
5	ar	0.8
6	trq	0.
7	trsS	0.8
8	trsA	0.2
9	qr	0.75
10	fps	0.073
11	fpq	0.086905
12	fns	fns
13	fnq	fnq
14	cs	5.
15	cf	100.
16	g
17	tMax	360
18	time2SlowTest	4
19	timeSlowTestWait	3

In[]:=

results=simul;

The chart below shows the state of each agent (represented by a column) through time (represented by each row). We can observe a significant proportion of the population getting infected.

In[]:=

results["MatrixPlot"]=mp[results]

Out[]=

S

,SQ

,E

,EQ

,IA

,IAQ

,IS

,ISQ

,R



Contents cannot be rendered at this time; please try again later

In[]:=

results["Plot"]=ListLinePlot[(results["TemporalData"]/@lAggCompartments)/VertexCount[results["Parameters"]["g"]],PlotLegendslAggCompartments,PlotLabelStyle["No Rapid Tests",Large],FrameLabel{"Time","Fraction Population"}]

Out[]=

Contents cannot be rendered at this time; please try again later

Tests performed were RT-PCR only.

In[]:=

ListLinePlot[100(Accumulate[results["TemporalData"]["SlowTest"]]/VertexCount[results["Parameters"]["g"]]),PlotLabelStyle["Cumulative Tests Performed per 100 agents",Large]]

Out[]=

Contents cannot be rendered at this time; please try again later

Results

Without rapid testing, tests are limited to people that have developed symptoms, and expanded to family members or other people with whom they might have been in contact with. This limited testing is not enough to eliminate the transmission within the population. Thus it requires societies to perform certain amount of lockdown or other restrictions to limit super-spreader events.

Experiment 2

In this case,the whole population is tested with the rapid test. We assume in this experiment that only 65% comply with the testing protocol, 35% of the population either refuse or avoid getting the rapid testing.

As in the previous experiment, the quarantine is not fully complied with. We have used a quarantine compliance rate of 75%.

In[]:=

aParameters["trq"]=.65;(*TestFractionofpopulation*)aParameters["trsS"]=0.;aParameters["trsA"]=0.;

In[]:=

results2=simul;

In[]:=

results2["MatrixPlot"]=mp[results2]

Out[]=

S

,SQ

,E

,EQ

,IA

,IAQ

,IS

,ISQ

,R



Contents cannot be rendered at this time; please try again later

In[]:=

results2["Plot"]=ListLinePlot[(results2["TemporalData"]/@lAggCompartments)/VertexCount[results2["Parameters"]["g"]],PlotLegendslAggCompartments,PlotLabelStyle["Rapid Protocol",Large],FrameLabel{"Time","Fraction Population"}]

Out[]=

Contents cannot be rendered at this time; please try again later

Epidemic is considered to have ended after 14 days without infected people.

In[]:=

endEpidemic=Last@First@SequencePosition[Differences[Accumulate[results2["TemporalData"]["Infected"]]["Values"]],ConstantArray[0,14]]

Out[]=

116

In[]:=

ListLinePlot[100((Accumulate[results2["TemporalData"][#]["Values"]][[;;endEpidemic]])&/@lTests)/VertexCount[results2["Parameters"]["g"]],PlotLegendslTests]

Out[]=

Contents cannot be rendered at this time; please try again later

In[]:=

results2//Keys

Out[]=

{States,TemporalData,Parameters,MatrixPlot,Plot}

In[]:=

results2["TemporalData"]//Keys

Out[]=

{S,SQ,E,EQ,IA,IAQ,IS,ISQ,R,FastTest,SlowTest,Susceptibles,Exposed,Symptomatic,Asymptomatic,Infected,Recovered}

Approximately 25-30% of the population will remain (or obey) the quarantine at a given time.

In[]:=

ListLinePlot[(Total[results2["TemporalData",#]&/@{"SQ","EQ","IAQ","ISQ"}]/1000)["Values"][[;;endEpidemic]],PlotLabel"Ratio of Population in Quarantine",PlotTheme"Detailed"]

Out[]=

Contents cannot be rendered at this time; please try again later

The majority of the quarantined agents are driven by the false positive rate results of the tests.

In[]:=

With[{pop={"SQ","EQ","IAQ","ISQ"}},BarChart[.001Transpose[(results2["TemporalData",#]["Values"][[;;endEpidemic]])&/@pop],ChartLayout"Stacked",PlotTheme"Detailed",ChartLegendspop/.rCompartments,PlotLabel"Fraction of population quarantine broken by compartment"]]

Out[]=

	Infected Symptomatic Quarantined
	Infected Asymptomatic Quanrantined
	Exposed Quarantined
	Susceptibles Quarantined

Results

Rapid, overactive testing + quarantine followed by individual that tested positive could allow to prevent the virus to continue to expand through the population while reducing the impact on the economy by quarantining those that have tested positive.

Conclusion

Nassim’s intuition of overactive testing that would be enabled by cheap fast and easy to administer rapid tests, seems to be a robust alternative to the current protocols being followed by the different public health officials at the county, state and country levels.

Initialization

In[]:=

Clear@simul

In[]:=

References

1 https://onlinelibrary.wiley.com/doi/full/10.1002/hed.26317

POSTED BY: Diego Zviovich

Answer

2 Replies

Sort By:

Moderation Team

Moderation Team, WOLFRAM

Posted 11 months ago

-- you have earned *Featured Contributor Badge* Your exceptional post has been selected for our editorial column *Staff Picks* http://wolfr.am/StaffPicks and Your Profile is now distinguished by a *Featured Contributor Badge* and is displayed on the Featured Contributor Board. Thank you!

POSTED BY: Moderation Team

Answer

Hamza Alsamraee

Hamza Alsamraee, Wolfram Research, Inc

Posted 11 months ago

I really enjoyed reading this Diego! I thought a recent post of mine might interest you, specifically in terms of schematic diagrams of compartmental models. Robert Nachbar actually made a CompartmentalModeling package that makes producing these types of diagrams very easy. More can be found here: https://community.wolfram.com/groups/-/m/t/2078640 I hope you find it useful!

POSTED BY: Hamza Alsamraee

Answer

Reply to this discussion

Reply Preview

Attachments

Remove Add a file to this post

Follow this discussion

or Discard

Group Abstract

Feedback