Abstract
Interventions to mitigate the spread of infectious diseases, while succeeding in their goal, have economic and social costs associated with them. These limit the duration and intensity of the interventions. We study a class of interventions which reduce the reproduction number and find the optimal strength of the intervention which minimises the attack rate of an immunity inducing infection. The intervention works by eliminating the overshoot part of an epidemic, and avoids a second-wave of infections. We extend the framework by considering a heterogeneous population and find that the optimal intervention can pose an ethical dilemma for decision and policy makers. This ethical dilemma is shown to be analogous to the trolley problem and we discuss how the dilemma can be avoided.
1 INTRODUCTION
Infectious disease epidemics have been suppressed and mitigated using a combination of non-pharmaceutical interventions (NPIs) such as lock downs, social distancing, mask wearing and contact tracing, and by immunising the population using vaccines. In the absence of vaccines, NPIs are the primary option. However, NPIs, and in particular, lockdowns, can have significant economic, mental health and social costs associated with them. Instead of protracted or repeated lock downs (as observed during the COVID-19 pandemic), a one-shot intervention has been suggested as a possible alternative. An intense but short-duration lockdown is imposed near the peak of the epidemic to stop the transmission during the overshoot phase of the epidemic and reduce the attack rate to the herd immunity threshold of the epidemic [1]. The overshoot phase is when the number of active infections start to decline (effective reproduction number is less than one), but a significant number of new infections are created. The overshoot is the difference between the attack rate and herd immunity threshold.
In this work, we explore an alternative strategy to achieve the same outcome through a prolonged but weaker intervention instead of a short and intense intervention. Such an intervention, if implemented early, will have the added benefit of reducing and delaying the peak of the epidemic as well, in contrast to the one-shot intervention [1]. If the economic and social costs of an intervention grow super-linearly with the strength of intervention, moderately intense but prolonged interventions may prove to be cost effective.
As with the one-shot intervention, the rationale of this strategy is to calibrate the intervention in such a manner that the attack rate of the mitigated epidemic is identical to the herd immunity threshold of the original epidemic. Therefore, when the intervention ends, there is no risk of further introductions developing into future epidemics. We show that this strategy is an optimal strategy for minimising the attack rate in the long term.
In the context of COVID-19 modelling, research on optimal interventions has attempted to include economic costs along with the objective of reducing infections: using detailed agent-based models [2] and fine-tuned intervention strategies [3], a balance is sought between socio-economic and health costs to minimise the total cost [4], or the claim that interventions reduce the economic well-being of a society has been challenged [5, 6]. An approach similar to the one presented here was suggested for COVID-19 in European countries [7]. Optimal interventions have also been studied as resource allocation problems where a limited stockpile of vaccine is available or a limited ‘amount’ of social distancing is acceptable and the objective is to find distribution of the intervention that minimises the reproduction number or a health-related objective function [8–10].
We do not include economic costs in an explicit manner in our model. The amount of reduction in ℛ0 can be interpreted as the cost – the higher the reduction in ℛ0, higher the social and economic cost of intervention. The calculations involved in finding the optimal strategy mainly rely on the knowledge of the basic reproduction number (or the next generation matrix).
We show that implementing an optimal intervention to minimise the attack rate could involve a moral/ethical dilemma for decision makers, which is analogous to the commonly known trolley problem [11, 12]. The dilemma arises as a result of transmission heterogeneity in the population.
In the following sections, we explain the modelling framework, results of our analysis, and conclude with a discussion of our modelling assumptions and the ethical dilemma that decision-makers could face.
II METHODS
We used deterministic SIR and SIR-like models to study the optimal intervention. In sub-section II A and II B, we explain the models used for a homogeneous population and for a heterogeneous population, respectively, in addition to describing the calculations for finding the optimal intervention. In II C we explain how an optimal intervention can be found if there is a delay in the start of the intervention.
A Homogeneous population
We use an SIR model with the variables s, i and r to represent the fractions of individuals in the total population who are susceptible, infected and recovered respectively [13, 14]. The size of the population is assumed to be constant such that, s + i + r = 1.
In this case, the attack rate of the epidemic is completely determined by the basic reproduction number ℛ0 and can be obtained using the following equation [13, 14]:
where s(t1) and s(t2) are the fractions of susceptible and r(t1) and r(t2) are the fractions of recovered individuals in the population at time instants, t1 and t2. Using the conditions i(t1) ≈ 0, r(t1) = 0 and i(t2) ≈ 0, which describe the population at the start and end of an epidemic, the attack rate relation can be obtained
An intervention that reduces transmission would affect the basic reproduction number as ℛ0 → ℛ0(1 − c) where 0 ≤ c ≤ 1. In the case of a homogeneous population, the herd immunity is achieved when the fraction of susceptible individuals in the population is less than 
. Therefore, we substitute s(t2) = 1/ℛ0 and s(t1) = 1 and solve for c. The optimal reduction in the basic reproduction number is
We verify this analytical result in section III A by simulating an epidemic where ℛ0 is changed to ℛ0(1 − c) in the early stage of the epidemic and the intervention is switched off once the active infections (i(t)) decline to a negligible number.
B Heterogeneous population
In a heterogeneous population, individuals may be further stratified into groups. To represent fraction of individuals in the total population who belong to a group k, we use the variables sk, ik and rk such that sk + ik + rk = nk, where nk is the proportion of the population who belong to group k and Σk nk = 1. Heterogeneity in transmission characteristics can affect the behaviour of epidemics in a significant manner. Epidemics in populations with different transmission structures but identical reproduction numbers can have widely different attack rates. An epidemic in a heterogeneous population can be described by the following SIR-like model, assuming identical duration of infection for all groups and measuring time in the units of the average infection duration,
The term Bkl is the average number of infectious contacts that an individual in group l causes in group k. The next generation matrix, G [15], can be constructed for this system with entries
The term Gkl is the expected number of infections caused in group k by an infected individual in group l. The dominant eigenvalue of G gives the reproduction number of the system [15]. The attack rates for this model are given by [16, 17]
It should be noted that we are defining the attack rate for the heterogeneous population as the number of infections in a group as a fraction of the total population.
The recipe for finding the optimal intervention is similar to the homogeneous case. Instead of a single reproduction number, we have the next generation matrix G. The herd immunity threshold is computed by finding the proportion of susceptibles of each group, sk, which would ensure that the reproduction number (top eigenvalue of G) is equal to one. In this case, there could be (infinitely) many such states that would ensure herd immunity. From all the possible herd immunity states, we select the one which maximises Σk sk (or equivalently, minimises the attack rate). The optimal intervention transforms B to B′ such that B′ results in the selected state as the attack rate.
A crucial point to note here is that unlike with the homogeneous population, it is possible for certain elements of B to increase in order to minimise the attack rate. In other words, the optimal intervention corresponds to an increase in transmission among certain groups or among pairs of groups. Thus, the change in reproduction number can not be a measure of of the economic or social cost in such cases. Nonetheless, this leads to some interesting results which are presented in the next section.
Instead of minimising the total number of infections, one can also consider a weighted objective function which is a weighted sum of the number of infections in each group. Such an objective function is useful when the disease affects different groups in different ways. The optimisation problem of finding the state of the population which minimises a general objective function and fulfils the herd immunity condition can be solved semi-analytically for the case of two groups and is presented in the Appendix.
C. Delayed intervention
In the above sections, we have assumed that the basic reproduction (or the next generation matrix) is a known entity and therefore an intervention is implemented right at the start of the epidemic. Calculating the strength of the optimal intervention requires knowledge of the reproduction number, the intervention would have to start after the epidemic has been established and enough observational data has been collected to calculate the reproduction number. While the basic principle would still hold, a delay could change the strength of the optimal intervention. To find the optimal strength for a delayed intervention, we use the attack rate relation with the final state s = 1/ℛ0, i = 0 and an initial arbitrary state sL, iL at a time instant tL when the intervention begins. We replace basic reproduction number in equation (1) with ℛ0(1 − c) and solve for
Using a numerical solution of the SIR equations, sL and iL can be found and the above equation can be solved for c. The equation (9) reduces to equation (3) when sL = 1, iL = 0, and c = 1 when sL = 1/ℛ0 and iL > 0. If 
, c > 1 which is biologically meaningless and reflects the fact that the population is already below the herd immunity threshold.
III. RESULTS
A. Homogeneous population
Simulation of the SIR model differential equations confirms our assertion in equation (3). As shown in FIG. 1, a ‘weak’ intervention reduces the attack rate but does not reduce the overshoot to zero. A strong intervention, on the other hand, reduces the attack rate as long as the intervention is in place but a resurgence occurs as soon as the intervention ends. The final health outcome under the strong intervention is worse than (or at least comparable to) the weak intervention, while incurring a higher social and economic cost during the intervention. The resurgence occurs because the small number of infections and sufficient number of susceptibles remaining in the population lead to new infections after the intervention is lifted. An intervention which is strong enough to minimise the final number of infections, while avoiding a resurgence, is the one whose attack rate (during the intervention) matches the herd immunity threshold of the unmitigated epidemic. This is the optimal intervention.
(A) Simulation of four types of intervention for an epidemic with ℛ0 = 1.5 in a homogeneous population: (i) No intervention – leads to largest attack rate, (ii) Weak intervention – reduces the attack rate, (iii) Strong intervention – reduces attack rate during the intervention, but leads to a resurgence in infections once the intervention is removed, (iv) A moderate intervention but optimal – attack rate is same as the herd immunity threshold. (B) The existence of a global minima for the attack rate shows that an optimal intervention strength exists. The resurgence of infections under a strong sub-optimal intervention is subject to certain assumptions which are discussed in section IV.
B. Heterogeneous population
Introducing heterogeneity in the model opens up a space of interventions that is not seen in the homogeneous case. In the homogeneous case, the herd immunity threshold is defined by a single point, but in the case of a structured population, the threshold is given by a collection of points. This can be seen by considering the following: the condition required for reaching herd immunity is that the typical infected individual must not infect more than one individual. In the homogeneous case, one can randomly chose a sufficient number of individuals and immunise them to ensure that the number of infectious contacts is less than one. If the population is structured, the typical infected individual must not infect more than one individual, on an average. As long as the average number of infectious contacts is less than one, herd immunity is achieved, irrespective of how the immunisation has been distributed among the various groups in the population. Thus there are infinitely many interventions that lead to herd immunity and prevent resurgence. Out of all these possibilities, we define the optimal intervention to be the one which minimises the attack rate.
When the population can be described using two sub-populations or groups, the optimal intervention belongs to one of the following types: the first group is fully infected, the second group is fully infected, or none of the groups are fully infected. This creates the possibility that under the optimal intervention, the number of infections in one of the groups is larger than what would have occured in the unmitigated epidemic, subject to the structure of the population. In FIG. 2 we show an example in which this occurs. This leads to an ethical dilemma wherein a certain group in the population incurs a higher cost (due to an increased number of infections) in order to minimise the cost for the whole population. Thus, the non linearity of the infectious disease dynamics, combined with population structure, lead to an ethical dilemma for policy/decision-makers which is analogous to the well-known trolley problem (see FIG. 2). The trolley problem involves a setup in which a train is going to hit a group of people who are lying on the tracks. The train can not be stopped, but a lever can be pulled to switch the train onto a different track on which only one person is lying. The dilemma that is posed by this situation is whether it is ethical to save a few lives at the cost of one?
Population structure and non linearity of infectious disease dynamics lead to an ethical dilemma. (A) A plot of the infections in the two groups that constitute the population. The corresponding cost function and next generation matrix are shown. The red curve shows the herd immunity threshold, the point marked with a cross shows the attack rate without an intervention, the point marked with a circle shows the attack rate if the optimal intervention is used. (B) A comparison of the interventions when the cost functions are different. The first plot shows that the optimal intervention leads to an increase in the number of infections in the first group. It is an example of the ethical dilemma of implementing the optimal intervention, which is explained further in (C). The second figure shows the plot for the cost function when the first group is given twice the weight as the second group, which means prevention of infections in the first group takes precedence over the second group. In this case, the intervention reduces the infections in both the groups. Thus, a biased cost-function may allow us to avoid the ethical dilemma. (C) This figure shows that the ethical dilemma involved in implementing optimal interventions is analogous to the well-known ‘trolley problem’. If the decision-maker does not act, the incoming epidemic (represented as a trolley) is going to cause many infections. If the decision-maker implements an optimal intervention (switches the tracks), the number of infections in the total population is minimised, but someone who otherwise would have been safe from infection, becomes infected (shown by the increased number of infections of group (II)).
In our case of optimal intervention, we are dealing with infections rather than deaths. So the dilemma maybe avoided by considering the cost of an infection for each group. Diseases often lead to a worse health outcome (mortality rate, hospitalisation rate, chance of leading to chronic conditions etc.) in certain groups of the population (the elderly age groups for instance). Instead of minimising the attack rate of the epidemic (which is the sum of attack rates in each group), it may be more prudent to minimise a cost function which is a linear combination of the attack rates in the groups, such that a group with a worse outcome of infection is given a higher weight in the cost function. For the example shown in FIG. 2, the ethical dilemma can be avoided using such a weighted cost function, although this is not always the case.
C. Delayed intervention
We calculate the optimal strength of the intervention and simulate the model to confirm the mathematical analysis in Section II C. Using equation (9) we observe that the strength of optimal intervention increases in a super-linear manner with the duration of delay. The results are presented in FIG. 3. As the population approaches the herd-immunity threshold, the strength of intervention approaches one – corresponding to the one-shot intervention [1].
Left: The time series of cumulative infections (i + r) for optimally controlled epidemics for a range of delays in the intervention. The vertical dash-dot lines show the time at which the intervention starts, and its corresponding time series is represented by the same color. Right: The strength of optimal intervention c plotted against the delay in intervention tL using equation (9). The super-linear increase in c shows the need for an early implementation. Parameters: Homogeneous SIR model with basic reproduction number ℛ0 = 1.5 and γ = 1. Intervention is implemented for a duration of 50 time units.
IV. DISCUSSION AND ONGOING WORK
In this work, we have examined a strategy of optimal interventions which uses the epidemic to cause just enough infections to induce herd immunity, eliminate the overshoot, and prevent future introductions from becoming epidemics. In addition to minimising the attack rate, this intervention would also slow down the growth of the epidemic and reduce the peak, which allows time to develop treatments and increase healthcare capacity.
In the case of heterogeneous transmission, our results indicate that the optimal strategy may require increasing the transmission in some of the groups and decreasing it in others, in order to minimise the attack rate for the whole population. This is analogous to the trolley problem, and it calls for a discussion around the ethics of subjecting certain groups to a higher rate of disease incidence, and the feasibility of this policy. If increasing transmission in certain groups is not viable either due to operational reasons or ethical considerations, herd immunity can still be achieved (and resurgence prevented) by reducing transmission in all groups. We have also presented some scenarios in which a weighted cost function avoids the dilemma. Our work shows that even without an explicit consideration of economic and social costs of an intervention, there are challenging questions to be answered for the first order problem of minimising the attack rate.
In addition to the ethical dilemma shown through our modelling here, interventions which require increasing transmission prompt an ethical discussion in relation to disadvantaged groups. Cultural, economic and social conditions factor into the contact structure of any human population – a high number of contacts due to living in close spaces, a high susceptibility to infection due to preexisting health conditions, or poor access to healthcare facilities, etc. Mathematical models of epidemics can throw light on possible choices of policy and may even help us pick the ones that lead to optimal outcomes. But the decisions made by policymakers are intertwined with political will, their popularity, and social attitudes. These eventually determine whether a particular intervention is favoured by a decision-making body [18, 19]. Disadvantaged groups, across the world, do not exercise sufficient political power to represent their interests in decision making bodies. In such a case, a decision-making body may find it convenient to subject a disadvantaged group to a higher attack rate in order to decrease the net attack rate for the whole population and achieve herd immunity, even if it does not minimise the net attack rate. The intervention strategy presented here, always carries such risks with it; and representation of disadvantaged groups thus becomes essential.
In addition to the ethical concerns with the optimal interventions presented here, there are some practical limitations to this strategy. There would be a natural tendency for individuals to protect themselves from getting infected even if interventions are not in place, so asking individuals to increase their transmission may not be a feasible strategy [20]. The optimal interventions could require a group of individuals to fully isolate themselves from the rest of the population. Such interventions are difficult to implement, as there would always be a small possibility for infections to be introduced into the isolated group [20]. If the transmission in other groups is increased, it would imply a larger chance of introduction into the isolated group.
We have assumed an SIR structure for disease progression in an individual. But, as long as the disease can be reasonably described by a model in which individuals do not become susceptible after getting infected, we would expect our results to be valid. A crucial detail that we have ignored is the stochastic and discrete nature of disease spread since it can capture the elimination behaviour of outbreaks, i.e., it can incorporate the difference between existence and absence of infections. The deterministic assumption and the use of continuous variables in our model means that after an intervention is over, the small number of infections present in the population will lead to another epidemic if herd immunity is not achieved (shown in FIG. 1, strong intervention). This however, is one of the possible outcomes. It is possible that the intervention completely eliminates all infections in the population, in which case a new epidemic does not result from any residual infections. However, even in this case the population remains vulnerable to an epidemic due to lack of herd immunity. Thus, a new epidemic can occur if new infectious individuals are introduced into the population. Another possibility is that the epidemic may get established with a delay due to the stochastic dynamics. Factors around contact-tracing and surveillance capacity (to eliminate the disease) and travel restrictions (to prevent introduction of new infections) are important for the selecting the optimal policy response, in addition to the results presented here.
The choice of two sub-populations was made to demonstrate the ethical dilemma arising out of the population structure, in the simplest manner. Generally, the population will be described by N sub-populations (or an N × N next generation matrix) with N increasing as the resolution of the information about the population increases. Therefore, describing an N -group structured population using an M -group model (where M < N) amounts to using a coarse-grained description of the system. We have shown that the optimal interventions in homogeneous population and in a structured population with two homogeneous sub-populations can be widely different. Thus, the optimal intervention can change if the amount of information about the structure changes. We intend to follow-up on this work by attempting to characterise how the optimal intervention changes as the scale of observation is changed, if the interventions could fail and if they can be robust.
Data Availability
All data produced in the present work are contained in the manuscript. Software used for analysis is available at https://github.com/Joel-Miller-Lab/optimal-intervention
VI. SOFTWARE REPOSITORY
The software used in this paper is available at https://github.com/Joel-Miller-Lab/optimal-intervention.
V. ACKNOWLEDGEMENTS
We thank Michael Meehan, Sarah Becirevic and Somya Mehra for helpful discussions. This work was supported by a La Trobe University Graduate Research Scholarship and a La Trobe University Full Fee Research Scholarship for Pratyush Kollepara, and startup funding from La Trobe University for Joel Miller.
APPENDIX FINDING THE OPTIMAL INTERVENTION
The generalised cost function depends on the weight vector a = (a1, a2, …) and the attack rate vector r = (r1, r2, …) and is given by:
The vector r is constrained such that ρ(G) = 1, where G is the next generation matrix and ρ(·) is the dominant eigenvalue.
We present a semi-analytical solution to this optimisation problem when there are two groups in the population. Minimising the cost function C is equivalent to maximising the function f :
while the constraint is given by equating the eigenvalue of the next generation matrix to one
where ∆ = det B. Values of s1 and s2 that satisfy this equation may correspond to cases where the non-dominant eigenvalue is unity. Therefore, the dominant eigenvalue needs to be computed numerically in order to impose this restriction.
To solve this optimisation problem we define a new function h = f + λg, where λ is a Lagrange multiplier. The partial derivatives of h with respect to s1, s2 and λ will be zero at points of interest and the three equations will be solved simutaneously to obtain
For a given matrix B, group sizes n1 and n2, and costs a1 and a2, the function f will be evaluated for the two points given by 
and the four boundary points. The six ‘types’ of possible extrema for the optimisation problem, are shown in TABLE I and in FIG. 4. Depending on the parameters, one of these points will maximise the function f subject to the constraint g = 0.
List of possible extremal points for a two-group population. The first column shows the solution type referred to in FIG. 4. The second and third column give the expression for the extremal values of s1 and s2. The quantities Bkl are the elements of the transmission matrix, {a1, a2} are the cost function weights. The quantity ∆ = B11B22 − B12B21.
A section of the parameter and optimal outcome space: The first column of each row shows the cost function, group sizes and the structure of the transmission matrix B. The next two columns are heat-maps which show the percentage change in the attack rate of Group 1 and 2 if the optimal intervention is implemented. The last column shows the type of solution obtained from the optimisation problem – the mathematical expression and description for the solutions can be looked up in TABLE I. The red coloring in the heat-map show regions of the parameter space where the optimal intervention would cause a worse outcome in one of the groups while minimising the cost function for the overall population. The model is identical for the first two rows, but the cost function is different. In the first row the cost function is the attack rate, while in the second row it is a weighted function of the two attack rates. The group sizes are 0.4 and 0.6 for all the rows.
The boundary extremal points represent scenarios where one of the groups is either fully infected or is fully protected from getting infected. Due to the large number (six) of parameters governing the optimal solution, we present graphical representations of the outcomes under the optimal intervention, in each group for sections of the parameter space in FIG. 4. We use a transmission matrix with a symmetric structure Bkl = bkbl((1 α)δkl + α) where α can be used to control the relative magnitudes of the diagonal and off-diagonal elements. The changes in the type of solution occur primarily based on the validity of the extremal points (
and 
) and which of the valid solutions minimise the cost. As an example, the straight line separation between type 5 (local minima) and type 2 (boundary point) and between type 5 and type 0 (boundary point), in the first row of FIG. 4, can be explained by imposing 
and 
on the expressions listed in TABLE I
