The optimal treatment of an infectious disease with two strains.

This paper explores the optimal treatment of an infectious disease in a Susceptible-Infected-Susceptible model, where there are two strains of the disease and one strain is more infectious than the other. The strains are perfectly distinguishable, instantly diagnosed and equally costly in terms of social welfare. Treatment is equally costly and effective for both strains. Eradication is not possible, and there is no superinfection. In this model, we characterise two types of fixed points: coexistence equilibria, where both strains prevail, and boundary equilibria, where one strain is asymptotically eradicated and the other prevails at a positive level. We derive regimes of feasibility that determine which equilibria are feasible for which parameter combinations. Numerically, we show that optimal policy exhibits switch points over time, and that the paths to coexistence equilibria exhibit spirals, suggesting that coexistence equilibria are never the end points of optimal paths.

Introduction

Several infectious diseases that pose a pressing problem in developing countries exist in multiple strains, including tuberculosis, HIV, Lyme disease, Human Papillomavirus, Helicobacter pylori and influenza (Balmer and Tanner 2011). Where resources are limited, a uniform treatment policy at the national level may be the most efficient way to tackle these diseases. This paper addresses the following question: what is the optimal treatment policy for an infectious disease that has two strains, which differ in their infectiousness? We take a standard model in epidemiology, the Susceptible-Infected-Susceptible (SIS) model, and extend it to include two strains of infection with different transmission probabilities. In the standard SIS model, individuals move between two states, susceptible and infected, based on exogenous probabilities. The probability of an individual catching a disease when he encounters an infected person depends on an infectivity or transmission parameter. This parameter is usually assumed to be homogeneous, which does not allow for policy differentiation if the disease exists in several strains.

The SIS model is often used to describe sexually transmitted diseases such as chlamydia. Chlamydia trachomatis is an organism that has three human strains: one strain causes trachoma, which is an infection of the eyes; a second strain causes sexually transmitted urethritis, often referred to simply as chlamydia, while a third strain causes the sexually transmitted lymphogranuloma venereum (LGV) (Byrne 2010). Chlamydia, caused by serovars D-K, is also the most common sexually transmitted disease in the United States and the United Kingdom, with over 1.2 million infections detected in the former and over 200 thousand in the latter.1 It can result in serious complications, such as pelvic inflammatory disease, which can lead to infertility if left untreated (Shaw et al. 2011). There are currently no vaccines for chlamydia, and treatment with antibiotics, usually with azithromycin or doxycycline, is necessary (Shaw et al. 2011). LGV is also treated with doxycycline (McLean et al. 2007), while trachoma is treated with azithromycin (Tabbara et al. 1996). Therefore, the optimal treatment mix for the various strains of chlamydia trachomatis is an important public health concern. It has been established in the literature that LGV strains are more infectious than trachoma strains, and that there may even be variation in virulence among trachoma strains (Pearce and Prain 1986; Kari et al. 2008). Reinfection with chlamydia is possible and superinfection does not appear to occur in practise (Heijne et al. 2011; Heiligenberg et al. 2014). Therefore, chlamydia trachomatis satisfies the assumptions of the model analysed in this paper, and as such, the results can provide guidance on the optimal treatment of the diseases that it causes.

We analyse the optimal treatment policy of a two-strain SIS disease, both analytically and numerically. We assume that the two strains differ in their infectivities, and neither superinfection nor eradication of the disease is possible. The policymaker can perfectly observe and target both strains with two separate treatment instruments. First, we provide theoretical results on the possible equilibria of this model. There are two sets of equilibria: non-boundary equilibria, where the two strains co-exist, and boundary equilibria, where one of the strains is endemic while the other is eliminated asymptotically. The feasibility of these equilibria depends on parameter conditions. Broadly speaking, when the infectivities of the two strains are very different, only the boundary equilibria can be achieved, and in particular those where the more infective strain prevails and the less infective strain is asymptotically eradicated. When the strains are quite similar, all equilibria are feasible. We demonstrate that in the coexistence equilibria, the more infective strain must receive more treatment than the less infective strain. These results lie in contrast to the equilibrium with no intervention, when the less infectious strain dies out.

We provide numerical examples where we simulate the model and look for optimal treatment policy along the path to the fixed points. We find that policy exhibits switch points along the path, implying that policy flexibility is beneficial, and that a static model would be too restrictive. We are unable to find a path towards a coexistence equilibrium that is optimal, in the sense of not being a spiral. We tentatively conclude that the coexistence equilibria can never be the end points of optimal paths, and that the boundary equilibria are always the optimal ones.

The ideas in this paper relate closely to the theory of competitive exclusion in ecology, which postulates that the coexistence of two or more species limited by one resource (in our case, two strains limited by one host population) is not possible (Armstrong and McGehee 1980). For example, Iggidr et al. (2006) demonstrate in the context of multiple-strain malaria that an equilibrium is either disease-free or has at most one surviving strain.

However, this is not always the case, provided some assumptions are relaxed, such as allowing for cycles in the persistence of species (Armstrong and McGehee 1980). In a model of renewable natural resources, coexistence of species is possible, with the stocks of the two resources higher in a mutualistic than a competitive equilibrium (Ströbele and Wacker 1991). In the case of an SIS model, coexistence of two strains is possible if one allows for superinfection and age structure (Li et al. 2010). Elbasha and Galvani (2005), in a model of two-strain Human Papillomavirus vaccines, show that when individuals are vaccinated for only one of the two strains, the equilibrium depends on whether the strains are synergistic (in which case vaccination also reduces the prevalence of the non-targeted strain) or antagonistic (where vaccination increases the prevalence of the other strain). Castillo-Chavez et al. (1999) demonstrate in a two-strain two-sex SIS model of a sexually transmitted disease that a unique coexistence equilibrium can occur under certain conditions. We show that in our two-strain SIS model, it is sometimes feasible for both strains to coexist, especially if they are not too different in their infectiousness. If one strain is much more infectious than the other, then it is only feasible to strive for the eradication of the less infectious strain, while allowing the more infectious strain to be endemic.

This paper builds on the economic epidemiology literature, which examines optimal intervention in tackling infectious diseases. Rowthorn et al. (2009) examine the optimal way to treat a disease with spatial dynamics, while Rowthorn (2006) considers the optimal treatment of a disease when policymakers have limited funds. There are important tradeoffs between treatment and vaccination as two instruments available to a policymaker (Rowthorn and Toxvaerd 2015). In general, a lack of intervention yields to inefficiently low uptake of vaccination (Chen and Toxvaerd 2014). While there is a rich literature that describes disease dynamics in the absence of intervention, the study of optimal intervention is still in its early stages, and this paper contributes to that literature.

The biology literature has shown the pervasiveness of diseases that exist in multiple strains, for which treatment and vaccination may have to vary. Lyme disease has at least 14 known strains.2 Seinost et al. (1999) present strong evidence that strains of Lyme disease have different infectivities. Other multiple-strain infections include Human Papillomavirus, which includes high-risk strains that are a leading cause of cervical cancer and low-risk strains that do not cause cancer, as well as tuberculosis, where some strains are particularly drug resistant (Castillo-Chavez and Feng 1997).3 Influenza is a multiple-strain disease: Truscott et al. (2012) argue that accurate modelling of influenza should not neglect the presence of several strains. HIV also exists in two forms: HIV-1 and HIV-2. Studies show that the less common HIV-2 strain is also less infectious than its counterpart for most of its infectious period.4 Accurate measurement of the transmission rate is the key to realistic simulations of HIV/AIDS prevalence rates for the USA and Africa (Oster 2005), which highlights both the importance of the transmission rate in describing an infection and the particular relevance of this paper to developing countries.

The paper proceeds as follows. Section 2 provides an overview of some theoretical concepts while Sect. 3 outlines the SIS model with two strains. Section 4 discusses the equilibria of the model and Sect. 5 provides examples of simulations. Section 6 concludes.

Relevant concepts

It is useful to review some relevant concepts before the central analysis of the paper. The SIS model belongs to the category of transmission system models, which are built on differential equations that describe the evolution of disease prevalence over time as a function of parameters. It has become standard to assume random mixing in models of this type: agents are equally likely to encounter any other agent in the population.5 There are two possible states: individuals are either susceptible or infected. They can move between the two states an unlimited number of times. Agents are homogeneous and the population is closed. If a susceptible agent encounters an infected agent, the susceptible agent becomes infected with a certain probability, termed the transmission probability. Typically, these models assume a homogeneous transmission parameter for all agents. The standard SIS model is in continuous time. At every time increment, a proportion of infected individuals is treated with a certain success rate of treatment (which can be interpreted as a rate of recovery). There is also the possibility of spontaneous or natural recovery.

Our model differs from the standard SIS model in two key ways: the disease has two strains with distinct transmission rates and the policymaker has two policy instruments. We make three further assumptions about the nature of the disease: transmission is through human-to-human contact rather than vector-borne, superinfection is not possible and agents can recover fully from the disease. SIS models are often used to describe sexually transmitted diseases (Keeling and Eames 2005), and chlamydia is an appropriate example of an SIS disease with multiple strains (Byrne 2010).

Optimal policy determines what proportion of individuals to treat at each point in time subject to the cost of treatment, the value of susceptible and infected individuals to social welfare, as well as the transmission dynamics. Studies focus on deriving equilibria and looking for policy switches over time. Policy generally takes the form of a Most Rapid Approach Path (MRAP)—the policy that takes the system the fastest way to the desired equilibrium. Another point of interest is Skiba points or curves, which delineate points where the policymaker is indifferent between policies.6

An analysis of optimal treatment in the one-strain SIS model can be found in Rowthorn (2006) and Goldman and Lightwood (2002). We summarise their results here as a brief indication of what one might expect from a two-strain model. Two types of equilibria of the model exist: the policy instrument is either at a boundary, where no one or everyone is treated, or policy is interior, with partial treatment of the population. It is never optimal to treat partially; it can be shown that optimal policy will always take on one of the two boundary values. Whether it is optimal to treat everyone or no one depends on the parameters at hand.

The SIS model with two strains

The model has the following characteristics. The two variants of the infection are High (H, high infectivity) and Low (L, low infectivity). The more infectious variant H has transmission rate while the less infectious variant L is characterised by transmission rate .7 The two policy instruments are , which measure the proportion of the infectees in each strain that are treated. The success rate of treatment is . We assume that the policymaker has perfect information about the prevalence of each strain in the population and can therefore target treatment of each strain perfectly.8

Individuals can catch either strain at the outset. When infected, they transmit the strain that they themselves are infected with. There is a possibility of exogenous recovery at rate . If individuals recover, they are again susceptible to either infection strain. The proportion of the total population infected with the High strain is . The proportion infected with the Low strain is . To ease calculations, the total population is normalised to size 1: , where S denotes the population of healthy individuals (susceptibles). All parameters are strictly positive. The policymaker places value p on healthy individuals, value 0 on infected individuals and discounts the future at rate . She also faces a constant marginal cost c of treatment.

In order to choose the optimal policy mix of and , the policymaker maximises the social welfare functionsubject to the equations of motion for the two infection types: Thus, the population of individuals infected with H increases at rate , which is a direct consequence of the assumption of random mixing, and decreases at rate . The population of agents infected with L changes in a symmetric way. We further assume thatThus, we assume that the population begins with a positive level of prevalence of each strain. Assumption (4) states that the infection probability in an encounter between a healthy individual and a sick individual is higher than the probability of a treated sick individual recovering. This ensures that neither variant of the disease can be eliminated, even asymptotically, by treating all infected people (see Sect. 4.4 for a proof of this). Thus, for all t.

In order to provide a benchmark for the optimal policy analysis, we discuss two simplifications of the model. First, we consider what happens when the two strains are independent, so that agents can be infected with both strains i.e. superinfection is possible, and infection with one strain has no impact on infection with the other strain. In this case, the disease dynamics are governed by the following equations:The equilibrium prevalence of the two strains isThus, prevalence depends only on the infectiousness of that strain and the extent to which it is treated. For the same level of treatment, the total number of infections (where a person infected with both strains counts as two infections) is higher than in the case when superinfection is not possible.

Second, we consider what happens to the infection when no intervention is possible; that is, at all points in time (but the infections are no longer independent). In this case, the two Eqs. (2) and (3) becomeFor each strain to be in equilibrium, we require . As long as , this can never hold with strictly positive and when there is no intervention. This implies that the two strains can never co-exist. Instead, there are two possibilities. Either the Low strain will die out and the High strain will prevail, or both strains will die out. The Low strain dies out because the probability of infection in each encounter of a healthy individual with one infected with the Low strain is lower than in the healthy individual’s encounter with one infected with the High strain. As a result, the High strain spreads more quickly and eventually causes the Low strain to disappear. The High strain will prevail if and die out if . Assumption (4) ensures that and so the High strain will prevail in equilibrium when there is no treatment, with .

To determine the optimal treatment policy in the general model with treatment, the policymaker solves the current value Hamiltonian:As the control variables and enter the Hamiltonian in a linear fashion, there is a bang-bang solution. The first-order conditions describing the optimal choice of treatment at a given point in time are where we drop time subscripts in these and many subsequent expressions for ease of notation. These conditions imply that if for example, then it is optimal to treat no one infected with the High strain. On the other hand, if , then the policymaker is indifferent between any level of treatment of the Low strain. These conditions can also be interpreted intuitively. Optimal treatment levels depend on the shadow price () of having an additional infected individual in the population. This shadow price will differ for the two strains. Considering the H type as an example, the optimal value of will depend on how compares to .9 The latter is the price of treating an individual divided by the probability that the treatment will be successful: it is the effective cost of treating the individual. It is constant for all infectees in the population. If the shadow price exceeds the effective cost of treatment, then all individuals are treated because the price of having them infected is greater than the price of treating them. On the other hand, if the shadow price is less than the effective cost, the price of treating an individual is greater than the price to social welfare of allowing this individual to remain infected and no one is treated. Finally, where the shadow price and effective cost are equal, the policymaker is indifferent between treating and not treating the individual and any level of treatment is optimal.

To complete the solution, the equations of motion for the two co-state variables are given below: These demonstrate how the shadow prices change with prevalence levels. The reproduction ratio describes the number of new infections arising from one infection (Anderson and May 1991). If it is greater than one, the infection persists; otherwise, it dies out. In this model, the reproduction ratios of the High and Low strains areAs a result of Assumption (4), these reproduction ratios are greater than one even if everyone infected with either strain is treated, so that neither strain can be eradicated (except in the case with no treatment, when the High strain will eradicate the Low strain).

Fixed points

In the analysis of dynamic models, of primary interest are the equilibria or fixed points of the model. In this section, we discuss the fixed points of the system and under which conditions they are attainable.

Analysis of non-boundary fixed points

We propose that there are two types of equilibria: non-boundary fixed points and asymptotic fixed points. They are defined as follows.

Definition 1

A non-boundary fixed point (FP) is a solution satisfying Eqs. (2), (3), (5), (6), (7) and (8), as well as .

Definition 2

An asymptotic fixed point (AFP) is a fixed point of the degenerate system in which one of equals zero and towards which there are paths with that converge asymptotically to this fixed point. We look for the treatment values in the vicinity of this fixed point that ensure this asymptotic convergence. On the boundary where one of the strains reaches zero prevalence (e.g. ), the co-state variable for that strain (e.g. ) will be infinite, while the co-state variable of the other strain (e.g. ) will be finite. In this sense, it is not a fixed point of the four variable system .

The key difference between AFPs and FPs is that some variables at an AFP are not constant: they move towards a constant but only reach it in the limit. The following Lemma reduces the set of feasible non-boundary fixed points to three:

Lemma 1

At any fixed point, the proportion of treated must be greater than the proportion of treated: .10 Therefore, only the following combinations of treatment at the fixed point are feasible: . This implies that there are three potential non-boundary fixed points.

Proof

At a non-boundary fixed point, it must be the case thatThis impliesEliminating yieldsand hence

Since it follows that Alternatively,At a non-boundary fixed point and From (2) it follows that the total prevalence of the disease is given bySince it follows from (7) thatEliminating yieldswhich can be written asThe above equation determines uniquely as a function of , and :At a non-boundary fixed point Subtracting (8) from (7), using the fact that and rearranging yieldsLemma 1 shows that at any non-boundary fixed point, so the above equation implies that

Following Wagener (2003), the parameter space can be split into different regimes that mandate which fixed points are feasible for which parameter constellations. We define the constant K, which will be used to create three regimes:The value of this constant is increasing in and but decreasing in and . Intuitively, this implies that K will be high when the difference between the two infectivities is high and when the success rate of treatment is low but the natural rate of recovery is high. We will show that there are three regimes of feasibility: and . We distinguish between two cases of fixed points, and .

Case 1:

We shall show that this case requires that and In this case because takes an interior value. Hence from (14), and From (9) and (14) it follows that this fixed point is characterised by the following values of the control variables and co-state variables:Thus, is decreasing in but increasing in . It is also increasing in because a more infectious strain requires more treatment but decreasing in , as this implies a lower overall spread of the disease. Since the parameters must satisfy the following condition:Economic logic requires that which yields the additional condition:Logic also requires that , which yields the conditionIn practice, this places constraints on c and p, which are otherwise free parameters. The value of can be determined from the equation for total prevalence:Total prevalence depends positively on and negatively on ; there is no dependence on as this enters through the choice of . The success rate of treatment, , does not affect total prevalence because no one infected with the Low strain is treated; the dependence of the High strain on enters, again, through the choice of .

Case 2:

Using (10), we find that in this case,For this solution to be feasible, we must have and henceTherefore, two subcases of this fixed point can be distinguished, depending on whether or

Case 2a Since , this implies that . In particular, is increasing in but decreasing in . There is interdependence between the two strains: the more infectious is the High strain, the less the Low strain is treated. However, the more infectious is the Low strain, the more it is treated. The optimal choice of also depends intuitively on and : the faster the rate of recovery, whether through treatment or naturally, the lower is the optimal treatment level. Since , it follows that and the value of is determined as follows:Thus, this fixed point is characterised by the following values of the control and co-state variables:Inserting the above values of and into (12) determines uniquely:The value of can then be determined from (11):which is increasing in but decreasing in and . This is intuitive: the more infectious is the H type, the higher the prevalence of the disease in equilibrium; on the other hand, the more successful the treatment or the higher the rate of natural recovery, the lower the total prevalence in equilibrium. A point of interest is that does not affect total prevalence; this is because the level of is fixed and chosen based on the value of . The dependence of on ensures that does not affect total prevalence. Note that we require and to be bounded by zero and one, so that an additional condition for the feasibility of this fixed point iswhich will place constraints on p and c.

Case 2b Equation (15) implies that in this case Hence andTo ensure that and we require thatThus, this fixed point is described as follows:There is one degree of freedom. The shadow price is indeterminate within the range shown. Eliminating and from (12) determines as a function of . The value of can then be determined from (11). There is a continuum of feasible fixed points with total prevalence

Summary of non-boundary fixed points

This summary uses multiple asterisks (*), (**), and (***) to denote distinct types of fixed point. Note that feasibility requires that and The above are the only feasible non-boundary fixed points. Note there are no feasible non-boundary fixed points if These findings suggest the possibility of a Skiba case when and . Inserting control variables in (11) shows that fixed point (**) always has higher total infection than fixed point (*). This is intuitive as in the former, both treatment levels are lower. For ease of notation, we label fixed point (*) as (low total prevalence), fixed point (**) as (high total prevalence) and fixed point (***) as (intermediate total prevalence).

If and there is a unique non-boundary feasible fixed point The control variables are and

If and there are two feasible non-boundary fixed points: with the control variables and ; and with control variables and

If the control variables satisfy In this case, there is a continuum of feasible fixed points which lie on the line

Analysis of asymptotic fixed points

The asymptotic fixed points are defined such that one strain approaches eradication asymptotically, while the other prevails at a positive level. As a result, the behaviour of the system in the neighbourhood of the fixed point can be approximated by the behaviour of a one-infection system. This is because the behaviour of the system with small is very similar to the behaviour when . Naturally, this holds for close to zero as well.

The behaviour of a one-infection system has been studied by Rowthorn (2006) and discussed in Sect. 2. Only boundary values for policy are optimal. This is because the co-state variable is a function of prevalence only in this type of problem. This implies that there must be a one-to-one mapping between prevalence and the co-state variable. However, the co-state variable changes with time. Therefore, a path cannot go back on itself in order to achieve two different values of the co-state variable for one level of prevalence. Interior values of policy require a path that goes back on itself (Rowthorn 2006). Therefore, interior values of policy cannot be optimal because they violate the assumption of the co-state variable being a function of prevalence.

To analyse asymptotic fixed points, we analyse a one-infection system and perturb it around this equilibrium point slightly. There are two possible categories of fixed points: those were tends asymptotically towards zero, and those where tends asymptotically towards zero.

boundary

First, suppose that The system then reduces to the following maximisation problem:subject toOptimal treatment is governed by the following bang-bang condition:Finally, the co-state variable evolves according to the following function:At a Hamiltonian fixed point, and :If the fixed point is interior then and we can solve the resulting equations to yield a unique interior fixed point There are also fixed points with and Let be the corresponding levels of infection. It has been shown in the literature (Rowthorn 2006) that it is never optimal to go to . Thus, the possible candidates for an optimum are and . Let

First, note that the point cannot be reached even asymptotically from the interior. Consider any point where :Recall that implies that . Hence,For sufficiently small values of x and y,  and for any value , the right hand side of this equation is positive. Hence, the prevalence of the High strain is increasing away from zero rather than decreasing towards zero, and cannot be reached even asymptotically from the interior.

What about the point Consider any point where :As before, , so thatFor sufficiently small values of x and y,  the right hand side is positive. In particular, for small x and y the right hand side is positive for all ifIf this condition is satisfied, there is no path that converges to from the interior. On the other hand, for , the asymptotic fixed point is attainable if . Thus, there is a feasible asymptotic fixed point with

Suppose next that . By symmetry, there are two potential categories of fixed points, with or . Let be the corresponding levels of infection. When are these feasible?

Let and Note that the point can always be reached asymptotically from the interior. Consider any point where . Suppose that and . In this case, the evolution of followsThe right-hand side of this equation is negative for small x and y when . Since an interior value of is suboptimal, this implies there is one final asymptotic fixed point, which is feasible for all K:What about with ? Consider any point where The evolution of is governed by the following equation:The right-hand side is negative for all values of Since interior values of can be ruled out by Rowthorn (2006), this asymptotic fixed points can be reached for all values of K by setting or :Thus, when , the only feasible candidates for asymptotic equilibria are those with .

Summary of asymptotic fixed points

If , there are three asymptotic fixed points:If there are no interior points that satisfy the Hamiltonian conditions. Moreover, the boundary points with cannot be reached. In this case, the only feasible candidates for an optimum are boundary points with In particular, there are three possible cases with the following control variables:The central conclusion when is that optimal policy involves asymptotic elimination of the low infectivity strain. The possibility of the Skiba case indicates that it is not enough simply to compare fixed points—the initial prevalence of the disease will affect which fixed point is optimal.

, with , and

with , and

with , and

with , and

with , and

Skiba case. The optimum path goes to either or depending on the starting point.

However, the following Proposition shows that asymptotic eradication of both strains of the disease is never possible:

Proposition 1

Both variants of the disease cannot be simultaneously eradicated even asymptotically in equilibrium ifwhich hold as a consequence of Assumption (4).

In Appendix B.

This proposition shows that even if both AFPs are feasible, asymptotic eradication of the disease is not possible and at least one strain always prevails.

Extensions to the model

In reality, strains do vary in factors such as drug resistance and the cost of treatment. In the present model, this could be modelled by varying the parameters and c by strain.

Suppose that the natural recovery rate of strain H is and that of strain L is . Similarly, the success rate of treatment of the High strain is and of the Low strain is . We can make some tentative claims about what would happen to the equilibria of the model. First, it should be noticed that coexistence of the two strains without intervention is now possible. The two strains evolve according to the following differential equations:When there is no intervention, the strains can co-exist ifIntuitively, this condition requires that the rate of infection and rate of recovery of the two strains are proportional. Thus, the High strain spreads its infection more quickly but also returns healthy individuals back into the population more quickly, so that the Low strain can survive. Next, we can analyse what would happen to the non-boundary fixed points. Lemma 1 states that only those fixed points where the High strain is treated more are feasible. In fact, if is higher than , and is higher than , this result could be reversed. A non-boundary fixed point requiresso that the relationship between the two treatment levels isThus, depending on the values of and , could be higher than . This makes intuitive sense, as higher values of and make the High strain comparatively easier to contain, so that this may outweigh the fact that a higher transmission probability makes the High strain more difficult to contain. Finally, we can consider the implications for the asymptotic fixed points. Here, also, the main result may be reversed: namely that the only feasible equilibrium is one where the Low strain is asymptotically eradicated. This is easy to see by rewriting Eqs. (19) and (21) allowing for differences in and by strain, and performing a similar exercise for asymptotic eradication of the Low strain. When is much higher than , and similarly is much higher than , then it may be that the only feasible equilibrium is asymptotic eradication of the High strain. This is because the low success rate of treatment of the Low strain, and its low natural rate of recovery, make it impossible to eradicate, even asymptotically.

Varying cost by strain would not affect prevalence or the fixed points that can be attained by the policymaker, but would affect the choice of equilibrium to converge to. This is ultimately a numerical question. It seems intuitive to assume that : it is more expensive to treat those infected with the High strain. The effect on the optimal equilibrium would depend how these values compare to the homogeneous cost parameter c. It could be that the overall cost goes down for a given level of treatment (such as if is slightly higher than c, but is much lower than c), in which case the policymaker would be more likely to pick an equilibrium where more people are treated, such as .

Feasibility

Summary of regimes of feasible equilibria

In this section we plot the feasible equilibria for each regime of K. Recall that if , the set of feasible equilibria is . If , the set of feasible equilibria is . If , the set of feasible equilibria is . As K increases, the set of feasible equilibria falls: there are fewer ways in which the disease can be contained. This is because a higher K occurs when treatment is ineffective or when the strains differ substantially in their infectivities, making it difficult to control both of them.

As K increases, the minimum attainable equilibrium prevalence of rises. This makes intuitive sense, as K is increasing in and in the difference between and . Where this difference is very large, and the only feasible equilibrium is the AFP where tends asymptotically towards zero and prevails at a positive level. In contrast, when is low, and one can attain equilibria such as , where the equilibrium prevalence of can be close to zero.

The set of feasible equilibria can be depicted graphically. Figure 1 shows the set of feasible equilibria when . The feasible set when is shown in Fig. 2. Figure 3 shows the set of feasible equilibria when .11 The graphs clearly illustrate the rising number of feasible equilibria as K falls.

Fig. 1

The set of feasible equilibria when . The full circles denote the two equilibria and

Fig. 2

The set of feasible equilibria when . The full circles show the equilibria and . The dashed line denotes the continuum of equilibria of type

Fig. 3

The set of feasible equilibria when . The full circles denote the equilibria (on the x-axis) and (on the y-axis). The empty circles denote the feasible equilibria and for this set of parameters, with the latter having higher total prevalence

An intuitive explanation and some examples

It is useful to introduce some parameter values in order to better understand the policy implications of these results. There are, in total, five individual fixed points: three non-boundary fixed points (one of which is in fact a continuum of feasible fixed points with the same total prevalence across the two strains) and two asymptotic fixed points. Feasibility depends on the value of . In order for the interior fixed points to exist, we require . If and , then two non-boundary fixed points exist, which suggest the existence of a Skiba line. A further consideration is whether the assumption is satisfied, which ensures that the disease cannot be eradicated, even asymptotically.

Generally, the existence of a non-boundary fixed point with positive equilibrium prevalence of the disease requires the following conditions, where (a) to (c) ensure the existence of the non-boundary fixed points and (d) ensures that at this non-boundary fixed point, the disease is not eradicated.We can interpret these conditions in an intuitive way. In order for both strains of the disease to co-exist in equilibrium, we require that their infectivities are not too different. This is because if the High strain is much more infectious, then individuals will be much more susceptible to it relative to the Low strain, so that the Low strain will eventually die out (asymptotically). On the other hand, if the infection probabilities are similar but low, relative to the recovery probabilities, then both strains will die out, because the rate of reinfection is not high enough to ensure the survival of the disease. The magnitude of the natural recovery rate relative to the success rate of treatment is important because a high natural recovery rate will likely lead to the elimination of both strains of the disease, while a low value of will imply that treatment is ineffective (and a particular level of treatment is required to ensure that the disease survives at a positive, constant level in equilibrium). Finally, the condition on the discount rate arises because economic logic requires in the non-boundary fixed point ; otherwise, the marginal social value of an additional infected person would be positive. Intuitively, not discounting the future too highly ensures that the impact of an additional infected person is felt for a longer period of time, hence ensuring that the marginal social impact of this infected person is negative at the fixed point.

The relative difference between the infectivities, , is not too high.

The natural recovery rate is not too high relative to the success rate of treatment: is not too high.

The future is not discounted too highly: is high enough.

The infectivity probabilities, and , are high relative to the recovery probabilities, and

Let us consider two examples that illustrate these points. First, suppose that , : the High strain is very infectious, relative to the Low strain. In order for , we require : the treatment needs to be very effective, relative to the natural rate of recovery. In addition, we require , in order for the disease to not be eradicated. One possible set of values for is (0.004, 0.04). Thus, in the example of a disease where the High strain is very infectious, the only way that both strains can co-exist in a non-boundary fixed point is when both the natural rate of recovery and the success rate of treatment are extremely low. Note that this example also requires a value of , which is reasonable, given that the economics literature usually assumes a discount rate of 0.95 or 0.99 (the latter is implied by an annual interest rate of 4%). If natural recovery is somewhat higher, e.g. if , then we will end up in the situation where the Low strain is asymptotically eradicated, and the High strain prevails in equilibrium. If either the natural recovery rate or the treatment success are high, e.g. if , even if , then both strains will be (asymptotically) eradicated.

Second, suppose that we have a very effective form of treatment with , but individuals naturally recover from the disease only rarely, so that . In this case, the disease will be eradicated, because the condition can never be satisfied. Suppose that instead we have a disease where and . Here, we require , so that infection is almost guaranteed when a susceptible individual meets an infected individual of either strain. One combination of values of that ensures coexistence of both strains in equilibrium is (0.99, 0.95). This set of parameters requires , which is easily satisfied.

The general lesson from these examples is that where treatment is very successful, the disease needs to be highly infectious in order to prevail. An interesting finding is that the rate of natural recovery relative to the success of treatment is also important: coexistence of both strains is more likely when natural recovery is low, but the treatment is very successful. Finally, coexistence of both strains is more likely when they are similar in their infectivities. This makes sense intuitively, as the bigger the difference, the more likely that the infectivity of the High strain will cause the Low strain to die out.

Feasible policies along the path to the fixed points

Another area where restrictions on feasibility may apply is the choice of treatment policy along the path to the non-boundary fixed points. Policies along the path towards the fixed points will always be boundary policies (i.e. 0 or 1), as these are Most Rapid Approach Paths. We take each of the boundary policies in turn and examine the feasibility conditions required for and to converge to their equilibrium values if this is the chosen policy. Different feasibility constraints will exist depending on whether the system begins below or above each equilibrium; below means beginning at an initial level , which requires and to increase to their equilibrium values, while beginning above implies an initial prevalence , which requires and to decrease to their equilibrium values. Letting denote the policy , , the feasibility restrictions are summarised in Table 1. Derivations can be found in Appendix C.

Table 1

Feasible policies along the path

	I˙_H,I˙_L>0	I˙_H,I˙_L<0
P^00	1-τβ_L>I_H+I_L	1-τβ_H<I_H+I_L
P^10	1-τβ_L>I_H+I_L	1-τ+αβ_H<I_H+I_L
P^11	1-τ+αβ_L>I_H+I_L	1-τ+αβ_H<I_H+I_L
P^01	1-τ+αβ_L>I_H+I_L	1-τβ_H<I_H+I_L

The table is interpreted as follows. If the system begins below , and need to increase. Not all policies can achieve this at all times: the range where the policy will allow and to increase is given by the top left entry in the table. There is an upper bound on ; if this condition is not satisfied, then setting will cause the system to move away from an equilibrium. The conditions in the table are more clearly shown in graphical form in Fig. 4.

Fig. 4

Feasible policies depicted for the case when , where the arrows denote the range where the policy is feasible. At , the policies are (0, 1) and (1, 1). The empty circles denote the feasible equilibria and for this set of parameters

The message of this figure is that when one starts below an equilibrium, feasible policies to reach the equilibrium usually involve treating no one or just one strain. This allows total prevalence to increase. On the other hand, when one starts above an equilibrium, feasible policies that ensure the equilibrium is reached treat at least one group, in order to ensure that prevalence falls. Thus, and can be reached from above with the policies and (1, 1), while they can be reached from below with the policies (1, 0) and (0, 0).

Summary of analytical results

Table 2 summarises the analytical results—the feasible fixed points, their associated optimal policies and for which values of K they are attainable.

Table 2

Summary of analytical results

Equilibrium	fH∗	fL∗	Feasible
F^high	(β_H-β_L)β_Lτα	0	K<1,K<δα
F^low	1	1-(β_H-β_L)β_Hτ+αα	K<1
F^int	1	0	K=1
ALhigh	1	0	K<1
AHhigh	0	∈{0,1}	All K
AHlow	1	1	All K

Simulations

In this section, we provide examples of optimal policy under various parameter combinations. When , the candidate equilibria are and . Simulations can tell us, under particular parameter values, which equilibrium it is optimal to go to, and what is the best way to get there. We explore the case when in the first part of this section. Next, we consider the situation when . In this case, there is a wide set of feasible equilibria, and it is not possible to simulate an optimal choice of equilibrium because the derivations of optimal policy are only correct in the neighbourhood of each equilibrium. Instead, we suppose that it is optimal to go to the non-boundary fixed point , and simulate optimal policy towards this point. We find that, in all the simulations that we estimate, optimal policy involves a spiral. We argue that this suggests that the interior fixed points are unlikely to ever be an optimal target for the policymaker.12

Optimal policy when

The values used for the parameters in the simulations are given in Table 3.

Table 3

Parameter values for K > 1

Parameter	Value
β_H	0.95
β_L	0.4
τ	0.15
α	0.2
p	1
δ	0.111

Some relevant concepts ought to be explained at this stage. First, we examine both optimal fixed and optimal variable policy. Optimal fixed policy compares the social welfare from choosing various policies at the beginning and not changing them until equilibrium. Optimal variable policy uses the Hamiltonian optimality conditions (5) and (6) to evaluate optimal policy at each time increment, using a fourth-order Runge-Kutta procedure. Paths with optimal variable policy are described as Hamiltonian paths because they satisfy these Hamiltonian optimality conditions. In the case of variable policy, we examine both forward and backward paths. This means that, starting from an initial prevalence, we simulate the system both towards and away from the fixed point. This is carried out by setting time increments equal to and respectively. As a result, it is possible to describe a longer optimal path than if we were to only simulate the system towards equilibrium.

Fixed policy

In this first set of simulations, we assume that policy is fixed along the entire path to equilibrium. This is a special case that can guide intuition, as policy may not always have the benefit of full flexibility. We assume that the policymaker chooses , and compare social welfare V when or . When the system begins, the policymaker chooses whether to treat everyone or no one infected with the High strain. The optimal policy is the choice of that gives the highest value of social welfare.13 The simulations set . Each time increment can be interpreted as one day, which implies that the results simulate an infection evolving over 90 days. The initial prevalence of the Low strain is set to 0.1.

Results are shown in Table 4. Optimal policy depends on costs. Where costs are low, it is optimal to treat everyone infected with the High strain. As costs rise, it becomes optimal to treat no one with the High strain. There is a unique optimal policy in regions I and III; in region II, optimal policy depends on the initial prevalence of the High strain. If the High strain is very prevalent, it is optimal to treat no one. If it is not, it is optimal to treat everyone. This strongly suggests that there exists a value of that is a Skiba point in Region II. By varying by small increments in region II and evaluating social welfare at each initial value for the two candidate policies, we find that a Skiba point exists when . At this point, the social welfare from selecting or is equal.14

Table 4

Optimal fixed policy as cost varies

Region	c	fH∗	fL∗
I (low costs)	c<0.2875	1	1
II (intermediate costs)	0.2875≤c≤0.3006	0 or 1, depending on IH0	1
III (high costs)	c>0.3006	0	1

Variable policy

In this set of simulations we allow policy to vary optimally over time. We look for strictly Hamiltonian paths; that is, we look for paths that are optimal because they satisfy the first-order conditions on and . For example, the policy ) is Hamiltonian at a particular point in time if at . In order to evaluate the Hamiltonian conditions, the co-state variables are needed at each point in time. These can be calculated using the facts thatThese partial derivatives can be approximated by perturbing the infection levels slightly. Thus, for initial infection levels and , the initial values of the co-state variables are for small . In these simulations we set . For a given initial infection level, we find the initial values of the co-state variables and using (25) and (26), and find the corresponding optimal initial policy. We then allow the system to evolve, checking for optimal policy at each increment of time by re-evaluating and using the equations for and that were derived in Sect. 3. We also simulate the paths backwards, starting at the initial point but moving away from equilibrium.15 Together, the forward and backward paths give a fairly complete picture of optimal policy over a period of length 2t, with in the middle. In this set of simulations we set , so that overall we examine the behaviour of the infection over 60 days. We assume costs are and set initial prevalence at .

We find that the system tends towards the equilibrium where , so that and This is consistent with the findings in Table 4, as falls into Region of costs. The evolution of the two strains is depicted in Fig. 5. This figure shows both the forward and backward paths. We find that prevalence of rises from a low initial point at first exponentially, and then roughly logarithmically until reaching the fixed point. Attainment of equilibrium takes approximately 15 days. The evolution of takes a slightly different path. Prevalence first increases, and then decreases towards zero. This is because as increases, there are fewer susceptibles that can become infected with the Low strain, so that the Low strain can be asymptotically eradicated.

Fig. 5

Evolution of prevalence along the path towards equilibrium , with and . Time >0 is the forward path, while time <0 is the backward path

When examining the sequence of optimal policy, we find that optimal policy exhibits switch points. This is depicted in Fig. 6. The treatment of the High strain switches from one to zero just before the start of the forward path, while the treatment of the Low strain switches from one to zero earlier.

Fig. 6

Optimal policy along the path towards , with and . Time is the forward path, while time is the backward path

Thus, optimal policy in this particular example has three phases:Thus, when policy is variable, it is always optimal to set in equilibrium when . However, policy does not begin at this level: it switches from one to zero when the marginal cost of an additional infected individual in the population is low, relative to the cost of treating him. In addition, the asymptotic eradication of the Low strain is possible in the closing part of the path even when no one in the Low strain is treated. Since this is feasible, it is also optimal, as it is cheaper than treating unnecessarily.

In this section we address the case when . In this regime, it is possible to have equilibria where both strains of the disease prevail. The parameter values are given in Table 5.

Table 5

Parameter values for K>1

Parameter	Value
β_H	0.8
β_L	0.4
τ	0.15
α	0.2
p	1
δ	0.222

We simulate optimal policy in the region of one of the non-boundary fixed points. In particular, we assume that it is optimal to go to the fixed point , and simulate optimal policy backwards, starting from that fixed point. The evolution of prevalence is shown in Fig. 7.

Fig. 7

Evolution of prevalence along the path towards equilibrium , with 0.4870 and . This path is simulated only backwards

This figure clearly shows a spiral. These simulations suggest that non-boundary fixed points are never the end points of optimal paths. In every case examined, perturbing the fixed point solution and running the path backwards generates a high frequency cycle. Reversing the path from its end point generates an integral which is exceedingly large. Figure 8 shows the evolution of optimal policy along this same path. The treatment of the Low strain switches between zero and one many times. Therefore, we conclude that even in the case when , it is optimal to move to one of the boundary fixed points. It appears that coexistence of the two strains is never optimal.

Fig. 8

Evolution of optimal policy along the path towards equilibrium , with 0.4870 and . This path is simulated only backwards

Chlamydia trachomatis

We began this paper with the example of chlamydia trachomatis, which is a disease with three strains that satisfies the behaviour of an SIS infection. As accurate modelling of the evolution of chlamydia trachomatis has been done in detail in the literature (e.g. Heijne et al. 2011; Sharomi and Gumel 2009; De Vries et al. 2006), we focus here on the distinction between the various strains, and the calculation of the cut-off parameter K, which is the novel contribution of this paper. In Table 6, we present the parameters that characterise chlamydia (serovars D-K) and LGV (serovars L1–L3) in the context of our SIS model with two strains.

Table 6

Parameters characterising chlamydia trachomatis

Parameter	Value	Source
β_H (LGV)	0.60	N/A
β_L (Chlamydia)	0.55	Althaus et al. (2012)
τ	0.18	Geisler et al. (2008)
α	0.97	Lau and Qureshi (2002)
K	0.02
Eradicate High strain?	Yes
Eradicate Low strain?	Yes

There are challenges to estimating the right parameter values for diseases such as chlamydia and LGV. Due to a lack of studies on LGV, we use the recovery rate and treatment success rate of chlamydia. The transmission parameter depends on the number of contacts between the infected partner and the non-infected partner, and we use the per-partnership parameter estimated by Althaus et al. (2012). We were unable to find an estimate of the transmission parameter for LGV, but we know that it is higher than that for chlamydia, so we set it marginally higher at 0.6. Next, the spontaneous recovery parameter is challenging to measure because of the ethical considerations of not treating individuals diagnosed with chlamydia (Geisler 2010). Geisler et al. (2008) screened individuals for chlamydia and then screened them again when they came for a treatment visit, for which the median was 13 days later, where they found spontaneous resolution in 18% of patients. Finally, the treatment success probability faces the challenge of non-compliance (Bachmann et al. 1999). We use the measured success rate of treatment with azithromycin in a meta-analysis of randomised control trials, which is found to be 97% (Lau and Qureshi 2002).

These parameters imply two observations about chlamydia trachomatis in our model: first, the value of K is less than one, so that the full range of fixed points should be attainable—those that involve coexistence of chlamydia and LGV, and those that involve asymptotic eradication of chlamydia and positive prevalence of LGV, although our simulations suggest that coexistence is not the end of an optimal path. Second, and more importantly, both diseases should be possible to eradicate, if everyone who is infected is treated. The question then is why the prevalence of chlamydia, for example, is around 6% among people younger than 25 in the UK (Macleod et al. 2005). The likely reason is that chlamydia is often an asymptomatic disease: 75-85% women are asymptomatic, while over 50% of men are asymptomatic (Risser et al. 2005). Thus, the treatment parameter is far from equal to one. In fact, a rough back-of-the envelope calculation, assuming values of and the fact that total prevalence in equilibrium is given by Eq. (3) implies a treatment parameter of .16 This suggests that only 35% of individuals infected with chlamydia at any given time are treated.

Conclusion

We have derived analytically the fixed points in an SIS model with two strains, where neither superinfection nor eradication of the strains are possible, and the policymaker has access to two separate treatment instruments. We have characterised two types of fixed points: non-boundary fixed points, where the two strains co-exist, and boundary fixed points, which involve the asymptotic eradication of one strain. The feasibility of these fixed points depends on parameter values: we have defined the parameter K, which delineates three regimes of feasibility. For , the only attainable equilibria are those where the Low strain is asymptotically eradicated and the High strain prevails. For , the coexistence of the two strains is a feasible equilibrium.

The simulations present two interesting results. First, when , optimal policy along the path to equilibrium exhibits switch points. This implies that a fixed policy is suboptimal. Further, in the neighbourhood of the equilibrium in our example, it is optimal to treat no one. This still allows the eradication of the Low strain. Second, when , we simulate optimal policy in the neighbourhood of one of the non-boundary fixed points. We find that optimal policy exhibits a spiral. This suggests that this equilibrium can never be the end point of an optimal path, and hence that coexistence, although feasible, is unlikely to be optimal. This suggests that for all values of K, the welfare-maximising policy is one that leads to one of the asymptotic equilibria.

Further research should consider extending this model to include protection via vaccination as another instrument available to the policymaker.

Table 7

Index of notation and acronyms

Notation	Definition
j=H,L	High, Low infectivity strain
I_j(t)	The prevalence of strain j in the population at time t
Ij0	The initial prevalence of strain j
f_j	The proportion of individuals infected with strain j that are treated
β_j	The probability that a healthy person encountering a strain j infectee will become infected
α	The success rate of treatment
τ	The natural rate of recovery from the infection
c	The marginal cost of treatment
p	The value of a healthy individual in social welfare
δ	The discount rate
λ_j	The multiplier in the Hamiltonian on the equation of motion for infection j
F^high	Fixed point where fH∗=(β_H-β_L)β_Lτα and fL∗=0
F^low	Fixed point where fH∗=1 and fL∗=1-(β_H-β_L)β_Hτ+αα
F^int	Continuum of fixed points where fH∗=1 and fL∗=0
ALhigh	The asymptotic fixed point where H is asymptotically eradicated, fH∗=1 and fL∗=0
AHhigh	The asymptotic fixed point where L is asymptotically eradicated, fH∗=0 and fL∗=∈{0,1}
AHlow	The asymptotic fixed point where L is asymptotically eradicated, fH∗=1 and fL∗=1
Z	The set of feasible equilibria
K	Equal to β_H-β_Lβ_Lτα, the parameter that characterises regimes of feasibility
P^ab	The treatment policy f_H=a,f_L=b
AFP	Asymptotic fixed point
FP	Fixed point
MRAP	Most Rapid Approach Path