Immune selection suppresses the emergence of drug resistance in malaria parasites but facilitates its spread

Although drug resistance in Plasmodium falciparum typically evolves in regions of low transmission, resistance spreads readily following introduction to regions with a heavier disease burden. This suggests that the origin and the spread of resistance are governed by different processes, and that high transmission intensity specifically impedes the origin. Factors associated with high transmission, such as highly immune hosts and competition within genetically diverse infections, are associated with suppression of resistant lineages within hosts. However, interactions between these factors have rarely been investigated and the specific relationship between adaptive immunity and selection for resistance has not been explored. Here, we developed a multiscale, agent-based model of Plasmodium parasites, hosts, and vectors to examine how host and parasite dynamics shape the evolution of resistance in populations with different transmission intensities. We found that selection for antigenic novelty (“immune selection”) suppressed the evolution of resistance in high transmission settings. We show that high levels of population immunity increased the strength of immune selection relative to selection for resistance. As a result, immune selection delayed the evolution of resistance in high transmission populations by allowing novel, sensitive lineages to remain in circulation at the expense of the spread of a resistant lineage. In contrast, in low transmission settings, we observed that resistant strains were able to sweep to high population prevalence without interference. Additionally, we found that the relationship between immune selection and resistance changed when resistance was widespread. Once resistance was common enough to be found on many antigenic backgrounds, immune selection stably maintained resistant parasites in the population by allowing them to proliferate, even in untreated hosts, when resistance was linked to a novel epitope. Our results suggest that immune selection plays a role in the global pattern of resistance evolution.

Z j (t) = S(t)(k i I(t) + k a (A j (t) + C(t))), where I(t) is innate immunity at time t, A j (t) is adaptive immunity against strain j, C(t) is cross reactive 10 immunity from another strain, and S(t) represents immune saturation.

11
Adaptive immunity: Adaptive immunity is calculated independently for each strain at each time 12 step, based the density of the strain and the duration of exposure to that strain in time steps. At time t, let 13 A j (t) be adaptive immunity to strain j, D j (t) be density of infected RBCs of strain j, E j (t) be the duration 14 of exposure to strain j, and A j (t+1) be adaptive immunity to strain j at time t + 1. If the host is currently 15 infected with strain j, the change in adaptive immunity can be described as follows: The first term describes the growth of the adaptive immune response, with the specific growth rate, g a , a 17 property of individual hosts. The second term is a decay term representing antigenic escape due to variant 18 switching. The constant ζ is the exposure duration at which the decay effect is half of its maximum, δ 19 governs the shape of the relationship between exposure duration and decay of immunity due to antigen 20 escape, and τ is the decay of immunity due to antigenic escape. As duration of exposure increases, fewer 21 novel antigenic variants remain and thus, the rate of switching slows and antigenic escape decreases 22 1 over time. These terms combined produce a pattern of dampened oscillations, in both immunity and 23 parasite populations. Adaptive immunity rises in response to expansion of the parasite population, 24 consequently causing a decline in parasite population. In response, density-dependent growth of the 25 immune response slows, and antigenic escape becomes the driving force. As adaptive immunity falls, the 26 parasite population expands again, thereby repeating the cycle. As exposure increases, the decay term 27 contributes less until finally adaptive immunity only rises over time, eventually eliminating the parasite 28 population.

29
If the host is no longer infected with the strain, it is simply a process of decay, described as: Adaptive immunity is constrained between 0.001 and 1.

31
Next, we calculate the contribution of cross-reactive immunity to the total adaptive immunity ex-32 perienced by strain j. Cross-reactivity is calculated from the non-self strain with the highest adaptive 33 immunity, given here as strain h. Let A h (t) be adaptive immunity against strain h at time t and C j (t) be the 34 cross-reactive immunity experienced by strain j at time t. The amount of cross-reactivity between strains 35 is given as χ and total cross reactive immunity is capped by χ max . The contribution of cross reactivity is 36 as follows: Innate immunity: Changes in innate immunity are strictly density dependent and are strain indepen-38 dent. Let I(t) be host innate immunity at time t, D(t) be total infection density, G be innate growth and H 39 be innate decay, so that innate immunity at time t+1 is calculated as follows: As calculated, I(t) may exceed one, and that value will be used to calculate I(t+1). However, the 42 contribution of I(t) to total immunity is capped at one.

43
Saturating immunity: As the density of parasites in an infection increases, the immune system 44 becomes overwhelmed, which is known as saturating immunity. As the infection density approaches the 45 saturation threshold, the immune system's total killing power decreases [2]. At maximum saturation, 46 we assume that the immune system's killing capacity is only 85% of its level of activation, a value given 47 by α. The relationship between saturation and density is given as follows, where S is the efficacy of 48 immune killing relative to its maximum value, D(t) is the total antigen density (here, assumed to include 49 merozoites, infected RBCs, and gametocytes) at time t, and the constant η determines the shape of the 50 relationship between density and saturation:

52
In each time step, RBCs die and are replaced with p RBCs, with total number not to exceed RBC carrying 53 capacity K. Uninfected RBCs are modeled as a pool, rather than individually, and so have no age. Let U 54 stand for the uninfected RBC count and M represent infected RBC count. The dynamics of uninfected 55 RBCs can be represented as follows:

56
In an uninfected host, background mortality removes a fraction m u of RBCs.
In an infected host, uninfected RBCs are subject to bystander killing, given as β, with additional 58 uninfected RBCs killed for each infected RBC, i.e, The multiplier, β has been estimated to be between 1 and 19, and may change over the course of an 60 given as V, with 1 representing a treated host and 0 an untreated host, and ω represents treatment efficacy 67 for sensitive parasites. The number of infected RBCs of strain j one time step ahead is then given by: We assume gametocytes share epitopes with parasites at other life stages, so gametocytes experience 69 the same immune process as well as daily background mortality. Gametocytes are not targeted by drug 70 treatment. We track the number of gametocytes of a given strain, L j , over time with: where γ is the daily fraction of gametocytes killed, B j (t) is the number of gametocytes of strain j maturing 72 at time t, and Z j (t) is immune killing for strain j.  can last as long as a primary infection because the short primary exposure leaves potential for antigenic 92 escape. Infections remain subpatent throughout. Transmission is unlikely in either case-even though 93 asexual replication is ongoing, high levels of immunity mean that most gametocytes will be killed during 94 the maturation delay, before they are transmissible.

95
Reinfection with a novel strain produces a brief, dense, symptomatic infection (S3B Fig). The infection 96 can grow rapidly due to the weaker adaptive immunity conferred by cross-reactivity, but the rapid 97 parasite expansion also triggers the innate immune response, curtailing the infection. These infections 98 are transmissible, albeit much less so than primary infections. Although adaptive immunity still kills 99 most gametocytes before they become infectious, enough are produced that some survive to maturity. infections. Finally, in these simulations, transmission intensity was high enough that most hosts were 113 repeatedly reinfected. This would weaken any impact of variation in duration of chronic infections.

114
Drug treatment rates 115 When the rate of drug treatment was increased from the default 30% rate, equilibrium prevalence (S6 Fig) 116 and rate of resistance evolution (S7 Fig) both increased, but the results were qualitatively similar to the 117 default treatment rate. The highest treatment rate that could be sustained without eradicating malaria in 118 any replicate simulation (out of 20) was 40%. Any replicates that were eradicated were excluded from 119 analysis. At treatment rates greater than or equal to 60%, the one strain and 300 vector condition achieved 120 a mean equilibrium resistance prevalence greater than our threshold for ubiquity (75%). For one strain 121 and 1200 vectors, the minimum treatment rate for ubiquitous resistance was 70%. Notably, even with a 122 70% treatment rate, the prevalence of resistance in the one strain conditions was lower and less stable 123 than the prevalence of resistance in the 30 strain conditions under the default treatment rate. Interestingly, 124 with a treatment rate of 20%, the overall pattern was similar to higher treatment rates in all conditions 125 except the 30 strain, 1200 vector condition. Resistance was substantially suppressed in this condition, 126 despite a effective treatment rate equivalent to the one strain conditions (S8 Fig), indicating that strength 127 of selection for resistance was not sufficient to overcome immune selection.

128
Because reinfection with a previously-exposed strain produced a low density infection, the majority 129 were more similar between strain conditions. This indicates that the lower equilibrium prevalence of 137 resistance in one strain conditions was not due to differences in effective treatment rate. resistance without recombination at four strain mutation rates, the default of 1×10 -5 , a higher rate that 148 was equal to the genomic mutation rate (2.5×10 -5 ), and two lower rates (5×10 -6 and 1×10 -6 ).

149
Recombination had no effect on T ubiq at any strain mutation rate ( it is not clear that this result can be generalized to natural populations.

164
At lower strain mutation rates, transmission intensity has little impact on T ubiq . In these conditions, 165 T ubiq was similar to T ubiq with 300 vectors at the default strain mutation rate (S9 Fig). This can be explained 166 by considering strain diversity patterns within the population (S11 Fig). Strain diversity and population 167 immunity were relatively low with lower strain mutation rates, and further, were similar between 168 transmission intensities. Within infected hosts, novel strain mutations were less frequent, decreasing the 169 opportunity for immune selection. As a result, high transmission did not delay the evolution of resistance. 170 At the higher strain mutation rate, the relationship between transmission intensities was qualitatively 171 similar to the default strain mutation rate, indicating a contribution from immune selection.