key: cord-318900-dovu6kha
authors: Pitschel, T.
title: SARS-Cov-2 proliferation: an analytical aggregate-level model
date: 2020-08-22
journal: nan
DOI: 10.1101/2020.08.20.20178301
sha: 
doc_id: 318900
cord_uid: dovu6kha

An intuitive mathematical model describing the virus proliferation is presented and its parameters estimated from time series of observed reported CoViD-19 cases in Germany. The model replicates the main essential characteristics of the proliferation in a stylized form, and thus can support the systematic reasoning about interventional measures (or their lifting) that were discussed during summer and which currently become relevant again in some countries. The model differs in form from elementary SIR models, but is contained in the general Kermack-McKendrick (1927) model. It is maintained that (compared to elementary SIR models) the model is more faithfully representing real proliferation at the instantaneous level, leading to overall more plausible association of model parameters to physical transmission and recovery parameters. The main policy- oriented results are that (1) mitigation measures imposed in March 2020 in Germany were absolutely necessary to avoid health care resource exhaustion, (2) fast response is key to containment in case of renewed outbreaks. A model generalization aiming to better represent the true infectiousness profile is stated.

Construction of the model has been motivated in course of the analysis of intensive-care capacity expenditure to be expected from sector-specific lifting of restrictions. This former analysis used a budget-oriented argument to arrive at an indicative estimate of the resource expenditure, but did not analyze dynamics. (Concretely it assumed a constant rate of new infections.) A reasonable question to be posed is: If a sector was allowed to reopen, what would the trajectory of infections actually look like, when surely it is not a linear increase? Further, can parameters of the local transmission behaviour be derived from the aggregate observed numbers?

In the present text, a model capturing the dynamics of the number of infections is developed towards answering these questions. It deliberately contains only few parameters and is in fact not designed to a specific stage of the virus proliferation. Though models for tracing the trajectory of infectious diseases exist, for example the intuitive SIR model [Ken56] which is formulated as a system of scalar differential equations, we believe that physically more realistic descriptions are possible, which, moreover, lead to increased accuracy of estimated parameters.

We assume a homogeneous set of individuals which act as unwitting agents in the proliferation. We assume that infection spreads probabilistically from the infected (and still contagious) individuals to any other individual of the set, wherein we assume that each individual is connected randomly to others, but such that all individuals approximately have an equal number of neighbours (=: A1). (In graph-theoretic terminology, the graph of contacts between individuals is a random undirected graph where each node has about the same edge degree d.) No other assumptions are imposed on the global topology of interconnections. Individuals who were once infected cannot be infected again (=: A2). Finally, an assumption here made is that contagiousness lasts only for a duration t c , i.e. an infected individual is contagious for the period [0, t c ] after its infection and then not at all afterwards (A3). This simplified characteristic is motivated by results on infectiousness found in epidemiological and clinical investigations: In [HLWea20] , infection incidence data of the Wuhan area is examined and combined with clinical data to derive an infectiousness profile which has most of its weight located at about 7 consecutive days around the symptom onset 1 . [WCG + 20] recorded viral RNA load data in sputum, throat swab and stool and report of nine patients viral peak loads of 2.35 · 10 9 copies per ml sputum, declining rapidly starting from the first day of presentation in almost all patients, decreasing to 10 5 copies per ml within about 10 to 16 days after symptom onset. (A level of below 10 5 copies per ml sputum, combined with no symptoms and past day 10, has been regarded as warranting discharge of the patient from clinical care with ensuing home isolation.) [TTL + 20] (Fig 2) report viral load in posterior oropharyngeal saliva samples decreasing monotonously to below 10 4 copies per ml in day 21 after symptom onset, for the majority of 20 non-intubated patients (out of n = 23). Changes in population size due to non-disease effects will be ignored, instead N will be considered constant; similarly the changes in proliferation characteristic due to disease-related reduction of the population will be deemed negligible. Assumptions A1 to A3 will be "baseline" assumptions throughout the text and substantial deviations from them will be discussed in the appendix only.

We aim for a numerical formulation of the aggregate evolution in which the randomness is averaged out. For this, let N be the number of agents, and let at t = 0 the number of infected agents x(t) be given as x 0 < N . Before t = 0, the number of infected agents shall be zero. To develop the model incrementally, lets momentarily assume that all infected agents are contagious infinitely long. In a unit time interval, all infected agents are deemed to infect each of respectively d other neighbours -stochastically independently -with probability p. The expected total number of virus receivers, per unit time interval, then is x(t) · d · p. But not all receivers get infected because some are already infected. The share of non-infected receivers among all agents is (1 − x(t)/N ); therefore the expected number of new infections in unit time is

Approximating the model evolution as continuous process even at small time intervals 1 Caution in the usage of numbers from pure incidence analysis is required: As consequence of the way the raw data is obtained in [HLWea20] , only infectiousness around the moment of symptom onset is in fact fully observed. This is because earlier transmission are likely usually not properly associated to the real primary case because the primary case does not show symptoms yet. Later transmissions are simply inhibited because the primary case is put into quarantine. An epidemiological analysis of incidence data alone therefore necessarily is insufficient to determine "pure" infectiousness. To emphasize the distinction between "pure"/"medical" infectiousness and infectiousness after taking into account the population's sociocharacteristics (household structures, current mitigation policies), it is worthwhile to call the density of the latter an "infection incidence profile".

(reasonable given the size of the numbers involved), one concludes, under assumption of infinitely enduring contagiousness, that x(t) followṡ

for t ≥ 0, with x((−∞, 0)) = 0 and x(0) = x 0 . Obviously the function x(t) is non-decreasing. For incorporating the finite duration contagiousness, one determines the number of contagious individuals as the difference of the accumulated number of infected at time t minus the accumulated number of infected prevailing at the earlier time t − t c , because that share of agents had the infection already for at least duration t c , and will cease to be infectious at t. Consequently, the expected total number of virus receivers is refined towards

The model with the finite duration contagiousness thus reads, in expectation,ẋ

with initial conditions as before. Both differential equations respectively have a unique solution. 2

The here presented model is not representable by the elementary SIR models that involve only instantaneous evaluations of the state variables on the right-hand side of the differential equation, as e.g. equation (2) in [Ken56] (see [ZML + 20] for a current example of its usage). It is therefore also necessarily different for example from [MB20] . The reason for this is a restriction imposed by such formulations, namely that the individual's transition from infection to recovery is modelled using a rate of transition proportional to the number of infected individuals, which corresponds to a stochastic recovery occurrence and yields an exponential decay characteristic on average. It is known however that, in reality, the SARS-CoV-2 shows a rather deterministic disease progression with regards to infectiousness in time, leading to end of the infectiousness after about two to three weeks after begin of infection, based on cell culture (see earlier citations). This clinically supported characteristic is properly represented in equation (2), but not in elementary SIR models. On the other hand the here presented model is conceptually contained in the original (i.e. general) compartmental model of Kermack and McKendrick [KM27] (which involves a formulation using integrals; see comments in [Bra17] also), for example by setting there ψ = 0. 3 This holds also for the refinement given in section B.

The advantage of the here given formulation is that it allows for a mathematically relatively simple description while still fully allowing accomodation of the infectiousness characteristic in generalized form. This simplicity gives some room to incorporate other, hithertho 2 Instead of considering only a temporally finite and uniform infectiousness, more detail can be incorporated into the differential equation using a convolution term, as shown in appendix B.

unconsidered, effects into the model and still retain a model complexity which is amenable to simulation for parameter identification.

3 Analysis of the model and exploratory simulation

For later simulation, it is helpful to make use of the scale invariances inherent in the above differential equations. If one denotes the equation (2) parametrized with d · p and t c and N and initial value x 0 as "ODE(dp, t c , N, x 0 )", then we have the following fact: If t → x(t) is a solution to ODE(dp, t c , N, x 0 ), then t → x(at) is a solution to ODE(a · dp, t c /a, N, x 0 ) for any a > 0. This means we can restrict analysis for example to t c = 1 and vary only d · p and x 0 . The other scale invariance is described by "x(·) solution of ODE(dp, t c , N,

Instead of a further analytical proceeding, the above equation's evolution was examined via computer simulation, for various parameter choices d · p and initial values. The purpose is first to explore the general (i.e. not real-data matched) behaviour of equation (2) (next subsection), then to fit the parameters to observed real data (section 4). Throughout it was used N = 1.0, t c = 1 and a (forward Euler) discretization step size of 0.01 (corresponding to resolution=100 in code).

The below discusses general features of the model and its behaviour under parameter variations. This is for demonstration only, and arguments on the proliferation phenomena should be taken as schematic. (Whether the phenomena occur in the real parametrization is to be discussed in section 4.) Fig 1a shows the evolution behaviour for some arbitrary but temporally constant parameter set. The most striking feature at this graph is that the number of infections asymptotically does not reach the total number N of agents. Rather, the limit is a value x(∞) < N which depends on the d · p and the initial value. For comparison, the evolution of the number of infections as would arise when observing eqn. (1) [with same d · p parameter] is depicted as grey dashed line; in it, the x(t) converges to N independent of the choice of d · p. (In subsequent text, this will be referred to as "bounded exponential growth".) The reason for including this curve here and in following graphs is that it can give a hint on trajectories of future viruses that may have a much more extended infectiousness interval. In fact, this curve would result if infected individuals remained infinitely long infectious and were not quarantined.

In simulations, the dependence of the limit x(∞) on d · p appeared to be generally overproportional (see Fig 1b) . This is well-known behaviour also in the instantaneous-state models. On the other hand, the dependence of the limit on x 0 was linear or sub-linear. In instantaneous-state models, the limit does not depend on the size of the initiating jump of Noteworthy is (in both cases) that even though the same number of exogenously infected was used as initially, the contagion effect is much smaller. The reason for this is that already about one fifth of the population had been infected (thus was immune in this model).

6 All rights reserved. No reuse allowed without permission.

perpetuity.

preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in

The copyright holder for this this version posted August 22, 2020. . https://doi.org/10.1101/2020.08.20.20178301 doi: medRxiv preprint

We use here the number of reported CoViD-19 cases (as aggregated by the Robert-Koch-Institut [1]) as a proxy for the number of infections in Germany. 4 We fit parameters for the interval until beginning of May 2020, assuming that the evolution proceeded within two different parameter regimes: first a d · p corresponding to no restrictions, then a d · p corresponding to the restrictions posed by contact disencouragement and store closure. (The observational interval used for parameter estimation does cover only a few days of the time of obligatory indoor face mask wearing.) We can derive parameters and based on them predict the trajectory of infections way forward. Because of the simplicity of the examined model, there is the risk of a high model error existing. Therefore, at the present state of this text, such estimation can only serve to determine reasonable bounds on the parameters of the model, rather than to give a reliable forecast of expect number of eventual infections.

Fitting of parameters is here conducted manually, focussing on moments in the time series that are indicative of parameter changes. At the beginning of April 2020, the number of weekly new CoViD-19 cases stood at about 40000 in Germany. If we regard the modelling time unit to correspond to a real duration of 2 weeks (implying that each individual newly infected is non-contagious two weeks after and onwards), then we have a new-infections rate of 80000 individuals per such time unit which corresponds to an increment of approximately ∆x = 0.001 per unit time after normalizing to N = 1.0. Identifying the moment which was one week after the initial wider lock-down in Germany (i.e. around 29th March) as moment t = 3 in the modelling, parameters consequently need to be fitted such thatẋ(3.0+) = 0.001 (green line). (The t = 3.0 also implies that the model assumes around 6 weeks of initial evolution under a low-restrictions scenario, which matches the timeline of the outbreak in Germany approximately.) Fig 3a shows the trajectory of the system evolution using initially d·p = 1.42 and switching to d·p = 0.7 afterwards. Fig 3b shows the evolution if no parameter switch (i.e. no intervention) had happened at t = 3.0.

Note: The matching is overly simplified for the interval t ∈ [0, 3.0], leading to an overestimated x(t), since for example x(3.0) = 0.0025 -corresponding to 200000 individuals-, while the actually reported number was around 52550. In reality, the d · p must have been larger than 1.42 at the beginning of the interval, but on the other hand closer to (but above) 1.0 in the second half of [0, 3.0].

Assuming an infected individual occupies an intensive-care bed with ventilator (ICU) for one to two weeks, the ICU capacity in Germany currently is about 12500 to 25000 ICU cases per week. This allows for a maximum of 87500 to 175000 reported infections per week (assuming share of cases needing intensive care around 14.28%), i.e. 175000 to 350000 reported infections per two weeks. This in turn corresponds to a normalized increment of 0.0021875 to 0.0043750 per time unit (a horizontal line somewhere in the upper half of the graphs in Fig 3) . 

It is necessary to remark that the conclusion drawn in connection with Fig 2e and 2f - i.e. that a second outbreak of similar magnitude as initially would not effect a substantial increase in the accumulated number of infected individuals -cannot be affirmed for the current scenario (in Germany and elsewhere), since that number is rather about 0.25% to 0.5% of total population currently, rather than the 1/5 prevailing in the demo scenario in Fig 2e and 2f at the onset of the second outbreak.

The challenge with lockdown measures for the current corona virus is the following: When imposing them, they will show effect only if the basic reproduction number is pushed below 1 sufficiently enough. Then, after the number of infected individuals has eventually dwindled, a lift of the lockdown is tempting -however even a slight increase of R 0 above one opens the way to renewed catastrophic infections increase. One therefore has a binary evolution characteristic; to control R 0 by policy such that a steady stream of just managable new infections is maintained is daunting, and likely impossible (in practice) if a policy requiring a constant set of restrictions is targetted. The natural answer, at least from a theoretical point of view, is to consider phases of lifted restrictions interleft with repeated adaptively switched phases of more stringent restrictions or more stringent enforcement of existing restrictions. The need for such strategy is not in principle altered by the local aspect of transmission, except that switched lockdowns only need to be local and thus do not affect the whole population.

Another point that needs to be mentioned is that the graphs suggest that a future virus having infectiousness lasting much longer than the about two to three weeks for SARS-Cov-2 and also being as highly infectious would pose serious challenges for containment, because of resource exhaustion in the mid-stages of the pandemic. 8 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in 

In this study a novel model for virus proliferation dynamics was developed and with it the SARS-Cov-2 outbreak in Germany retraced on an aggregate level, using CoViD-19 case count data by the Robert-Koch Institute in Berlin. Elementary properties of the model were identified. Predictions by the model for different levels of mitigation measures were hinted at or stated in approximate manner, and put into context of available health care resources in Germany.

Future policy oriented work would need to address better understanding of fine-grained and adaptively activated mitigation measures, for which a spacial model should be favoured over purely aggregate models as the present one. Further, for purpose of improving parameter and state estimates, the issue of underreporting (i.e. #actual > #reported cases) must be taken into account appropriately. Ideally, one can develop an estimate for the factor of underreporting from more exact spacial analyses.

On the mathematical side, a more rigorous formulation of the instantaneous proliferation dynamics is desirable, which allows to link parameters of the aggregate model to well-defined elementary parameters and results in more systematic parameter estimation. The ultimate goal is to be able to estimate more local structure from the observed time series.

10 All rights reserved. No reuse allowed without permission. perpetuity. preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in The copyright holder for this this version posted August 22, 2020. . https://doi.org/10.1101/2020.08.20.20178301 doi: medRxiv preprint A Additional graphs A.1 Data series on daily newly reported CoViD-19 cases Figure 5 : A "smoothed" derivate of numbers of daily newly reported CoViD-19 cases in Germany published by [Rob] . The blue squares and green triangles series show (for comparison) the sum of daily new cases over a moving 7-day window. Orange diamonds and triangles show daily new cases after scaled with a weekday-specific weight factor to remove the weekly pattern seen in the original data. The weight factors were estimated from data corresponding to the squares and diamonds series, i.e. from the interval from 1st April until 6th May 2020. Germany imposed face-mask wearing in stores starting from 27th April and allowed certain (moderate) shop reopening starting from 4th May 2020. The "bend" at around 14th April is remarkable because no changes in measures were effected at that time or within the preceding one week.

B Refinement of the infectiousness mechanism, including a model generalization So far, a crude specification of the infectiousness has been used, putting focus on the main infectiousness interval of a few days. An additional aspect in the virus transmission which should be accounted for in a refinement is the transmission from longer lived remnants of the virus in otherwise cured individuals. For this, we imagine that individuals infected at time t 0 remain contagious until t 0 + t c2 with reduced probability, additionally to the previously used interval [0, t c ]. Concretely, let p 2 be the probability that an individual which has been infected for a duration exceeding t c but not exceeding t c2 , will transmit the virus in a unit time step. Withp 2 := p 2 /p the adjusted model equation then readṡ

x(t) = d · p · (x(t) − x(t − t c )) +p 2 (x(t − t c ) − x(t − t c2 )) · (1 − x(t)/N ),

11 All rights reserved. No reuse allowed without permission.

perpetuity.

preprint (which was not certified by peer review) is the author/funder, who has granted medRxiv a license to display the preprint in

The copyright holder for this this version posted August 22, 2020. .

Mathematical epidemiology: Past, present, and future

Temporal dynamics in viral shedding and transmissibility of COVID-19

Deterministic and stochastic epidemics in closed populations

A contribution to the mathematical theory of epidemics

Effective containment explains subexponential growth in recent confirmed COVID-19 cases in China

COVID-19: Fallzahlen in Deutschland und weltweit

Temporal profiles of viral load in posterior oropharyngeal saliva samples and serum antibody responses during infection by SARS-CoV-2: an observational cohort study

Katrin Zwirglmaier, Christian Drosten, and Clemens Wendtner. Clinical presentation and virological assessment of hospitalized cases of coronavirus disease 2019 in a travel-associated transmission cluster

Early prediction of the 2019 novel coronavirus outbreak in the mainland China based on simple mathematical model

since those individuals must be added to the instantaneous reservoir from which infections are generated. The equation is better written aṡIf we denote by i(t) the infectiousness profile, which shall describe the relative infectiousness of an infected individual 5 at time increment +t after the infection moment (relative to infectiousness at t = 0), then the above used specification for SARS-Cov-2 is expressed asIts derivative is (with Dirac notation) i = δ 0 − (1 −p 2 ) · δ tc −p 2 · δ t c2 . We therefore find that the model equation (6) in fact is generally best written aṡwhere " * " denotes the function convolution. A Dirac notation-free representation derivesHere the last integral signifies the well-known Stieltjes integral.Note: The infection from contaminated surfaces of objects can be represented in the same framework. This is because initially and during the evolution of the spread, viruses are on surfaces mostly there where infected individuals previously had been.