(C) PLOS One
This story was originally published by PLOS One and is unaltered.
 . . . . . . . . . .



Evolutionary safety of lethal mutagenesis driven by antiviral treatment [1]

['Gabriela Lobinska', 'Department Of Molecular Genetics', 'Weizmann Institute Of Science', 'Rehovot', 'Yitzhak Pilpel', 'Martin A. Nowak', 'Department Of Mathematics', 'Department Of Organismic', 'Evolutionary Biology', 'Harvard University']

Date: 2023-08 

Nucleoside analogs are a major class of antiviral drugs. Some act by increasing the viral mutation rate causing lethal mutagenesis of the virus. Their mutagenic capacity, however, may lead to an evolutionary safety concern. We define evolutionary safety as a probabilistic assurance that the treatment will not generate an increased number of mutants. We develop a mathematical framework to estimate the total mutant load produced with and without mutagenic treatment. We predict rates of appearance of such virus mutants as a function of the timing of treatment and the immune competence of patients, employing realistic assumptions about the vulnerability of the viral genome and its potential to generate viable mutants. We focus on the case study of Molnupiravir, which is an FDA-approved treatment against Coronavirus Disease-2019 (COVID-19). We estimate that Molnupiravir is narrowly evolutionarily safe, subject to the current estimate of parameters. Evolutionary safety can be improved by restricting treatment with this drug to individuals with a low immunological clearance rate and, in future, by designing treatments that lead to a greater increase in mutation rate. We report a simple mathematical rule to determine the fold increase in mutation rate required to obtain evolutionary safety that is also applicable to other pathogen-treatment combinations.

Funding: This work was supported by the Kimmel Foundation (YP) https://www.weizmann.ac.il/acadaff/awards-and-honors/institute-honors-and-prizes/Kimmel ; and the Minerva Center on Live Emulation of Evolution in the Lab (YP) https://www.minerva.mpg.de/centers/list/minerva-center-on-live-emulation-of-evolution-in-the-lab . YP holds the Ben May Professorial Chair. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. The authors did not receive any salary from any of the funders.

Copyright: © 2023 Lobinska et al. This is an open access article distributed under the terms of the Creative Commons Attribution License , which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Introduction

Nucleoside analogs are molecules similar in shape to naturally occurring nucleosides used by living organisms and viruses for nucleic acid synthesis. They are therefore readily incorporated into nascent DNA or RNA chains by viral polymerases. Many nucleoside analogs differ from natural nucleosides in key aspects which usually prevents further viral genome chain elongation. Some nucleoside analogs lack a 3’OH group which makes the viral polymerase unable to attach the next nucleoside to the growing chain. Others, such as Lamivudine, cause steric hindrance upon incorporation into the DNA or RNA chain [1–3].

Other nucleoside analogs do not prevent viral RNA elongation. Instead, they have the capacity to ambiguously base pair with several nucleosides, causing erroneous incorporation of nucleosides during the replication process. Thereby such drugs increase the virus mutation rate up to the point of lethal mutagenesis, a mechanism with foundations in quasispecies theory. This theory describes populations of replicating genomes under mutation and selection [4–9]. In addition to increasing the probability of lethal mutations, mutagenic treatment can also decrease the number of viable virions through lethal defection [10]. Lethal defection occurs when functional proteins synthetized from viable viral genomes are consumed for the packaging of defective ones that coexist in same cell.

Repurposing mutagenic antiviral drugs to treat Coronavirus Disease-2019 (COVID-19) has been suggested early on in the pandemic [11]. Molnupiravir, a prime example, seems to act exclusively through mutagenesis. Its incorporation into nascent RNA genomes by the viral polymerase does not result in chain termination, in fact, the viral RNA polymerase has been shown to successfully elongate RNA chains after the incorporation of Molnupiravir [12–14]. Molnupiravir switches between 2 tautomeric forms: one is structurally similar to a cytosine, the other is structurally similar to a uracil. Hence, Molnupiravir can base pair, depending on its form, either with guanosine or with adenosine [12,13]. Severe Acute Respiratory Syndrome Coronavirus 2 (SARS‑CoV‑2) is a positive-sense single-stranded RNA virus and its RNA replication proceeds in 2 steps. First, the negative-sense RNA is polymerized based on the plus strand, and the negative strand then serves as a template to synthetize positive-sense RNA molecules [15]. Hence, the incorporation of Molnupiravir during the first step of RNA synthesis gives rise to an ambiguous template: positions where Molnupiravir was incorporated can be read by the RNA-dependent RNA polymerase as either guanosine or adenosine. This causes mutations in the progeny RNA compared with the parental RNA, possibly up to the point of the “error threshold” and death of the virus [12–14], see Fig 1A. For a discussion of error threshold and lethal mutagenesis, see [6,7,16–20].

PPT PowerPoint slide

PNG larger image

TIFF original image Download: Fig 1. (A) Mechanism of action of molnupiravir. SARS‑CoV‑2 has a positive-sense single-stranded RNA genome, represented schematically in (1). Its replication proceeds in 2 steps: first, the synthesis of a negative-sense template strand (2), which is then used to synthesize a positive-sense progeny genome (3). Molnupiravir (M) is incorporated against of A or G during the synthesis of the negative-sense template strand (2). When the template strand is replicated, M can be base-paired with either G or A. Hence, all A and G in the parent genome become ambiguous and can appear as A or G in the newly synthetized positive-strand genome. C and T are not affected by molnupiravir during the synthesis of the template strand, (1) to (2), but can be substituted to M during the synthesis of the progeny genome from the template strand, (2) to (3). M can then base-pair with A or G when used as a template; see (3) to (4), which can cause A->U and U->A transitions in the final progeny genome (5). (B) Virus dynamics within an infected person. Wild type (x) and mutant (y) replicate at rate b and quality q = 1−u. The per base mutation rate, u, is increased by treatment with molnupiravir. Both wild type and mutant need to maintain m positions to remain viable. Mutating any of n positions in the wild type results in a mutant. In the beginning of the infection, the adaptive immune response is weak, and virus is cleared at a rate a 0 which is less than b. After some time, T, the adaptive immunity is strong, and virus is cleared at the rate a 1 which is greater than b. (C) Graphical summary of the influence of mutagenic drugs on virus mutants. White circles represent wild type, beige circles viable mutant, and black circles dead virus. When the mutation rate is low, few viable mutants and few lethal mutants are produced. Most mutations occur when the virus load is already high; hence, they have little influence on subsequent generations. For intermediate mutation rate, the total virus load declines but the amount of viable mutant increases. When the mutation rate is high, both the virus load and the amount of viable mutant decline. SARS‑CoV‑2, Severe Acute Respiratory Syndrome Coronavirus 2. https://doi.org/10.1371/journal.pbio.3002214.g001

As noted before, the intended antiviral activity of Molnupiravir resides in its capacity to induce mutagenesis and hence reduce virus load. Yet, this very property which confers to Molnupiravir its desired antiviral effect might also enhance the capacity of the virus to develop drug resistance, immune evasion, infectivity, infectiousness, or other undesired phenotypes. Thus, a mathematical analysis should weigh the desired and potentially deleterious effects of mutagenesis drugs in general and of the present virus and drug in particular.

Mathematical theory has established the concept of an error threshold, which is the maximum mutation rate compatible with adaptation (or survival) of a population of replicating agents [9,21,22]. But a theoretical analysis is still missing to evaluate the risk of emergence of variants of concern (VoC) or simply viable mutants following mutagenic treatment. In the context of the COVID-19 pandemic, the mutagenic potential of Molnupiravir leads to concerns about accelerating SARS‑CoV‑2 evolution. VoC can include mutants that are resistant against vaccination or antiviral treatments as well as mutants with enhanced transmissibility or lethality.

Since the beginning of the COVID-19 pandemic, mathematical modeling has been used to establish vaccination strategies minimizing the risk of emergence of resistant mutants [23–27], predict the epidemiological spread of SARS‑CoV‑2 [28–32], and optimize the guidelines concerning isolation of contact and positive cases [33,34]. The emergence of resistance against treatments has been investigated for other viruses, including human immunodeficiency virus (HIV) [35,36], influenza A [37–40], and hepatitis C [41].

Molnupiravir has been found effective in inhibiting the replication of SARS‑CoV‑2 in ferrets, mice, and cultured human cells [42,43]. Following these promising results, Phase 2 and then Phase 3 clinical trials were conducted and concluded that Molnupiravir is safe and reduces the risk of hospitalization or death by about 50% [44–47]. The recommended dosage is 800 mg twice daily, for 5 days, and within 5 days of symptom onset [48]. It is especially recommended for individuals at high risk for disease progression to severe symptoms and death.

While drug’s physiological safety is a cornerstone of pharmacology, we explore here a new aspect of drug safety. We define a treatment as “evolutionarily safe” if its use does not increase the rate of generation of mutants in a treated patient beyond the risk expected in an untreated patient. This notion is especially relevant for drugs that work through lethal mutagenesis.

In this paper, we analyze the case study of the increase of the evolutionary potential of a virus (here: SARS‑CoV‑2) under mutagenic treatment (here: Molnupiravir treatment). In particular, we ask if the wanted effect of limitation of virus load by the drug could be accompanied by an unwanted enhancement in the rate of appearance of viable mutants or VoC due to increased mutagenesis. We construct a mathematical framework describing the increase and decrease of the virus load after infection and derive expressions for the total amount of wild type and mutant produced by individuals during the course of an infection. We use empirical data on COVID-19 and bioinformatic data on SARS‑CoV‑2 to estimate key parameters, including infection progression within the body amidst response of the immune system and the number of potentially lethal positions in the genome.

We find that the Molnupiravir-SARS‑CoV‑2 treatment is situated in a region of the parameter space that is estimated to be narrowly evolutionarily safe. Evolutionary safety is expected to increase in patients with decreased immunological viral clearance rate. Our analysis also shows that evolutionary safety increases with the number of positions in viral genome positions that are lethal when mutated. Crucially, evolutionary safety could be improved by obtaining higher increases in the mutation rate under treatment that provides a clear direction for future drug improvement. We derive a simple mathematical formula that determines the evolutionary safety of a drug given the pathogen’s mutation rate with and without treatment and the number of positions in the pathogen’s genome that are lethal when mutated.

Description of the model After infection with SARS‑CoV‑2, virus load increases exponentially until it reaches a peak after a median of about 5 days [49]. During this growth phase, the action of the immune system is insufficient to counterbalance viral replication. Subsequently, the immune response gains momentum and infection enters a clearance phase. Now virus load decreases exponentially until the virus becomes eliminated about 10 to 30 days after initial infection [49,50]. In some immunocompromised individuals, viral clearance can take many weeks [51,52]. However, some argue that the isolation of infectious virus is rare after 20 days postinfection [53]. The values for the time to peak of the virus load, time to clearance, and the virus load at peak from several sources are summarized in Table A1 in S1 Text. We provide an overview of related published literature, in the S1 Text (see section “Relationship to previous literature”). In our mathematical formalism, we describe the evolution of a virus within the body of a single human host by following the abundance of 2 viral types: wild type, x, and mutants, y. Both x and y replicate with birth rate b and replication quality q = 1−u, where u is the mutation rate per base. The mutation rate can be altered by the administration of a mutagenic drug. The virus genome contains m positions, all of which must be maintained without mutations in order to generate viable progeny. In addition to those, we consider n positions, such that even a single mutation in one of them gives rise to a mutant virion, y. In this paper, a mutant virion can be any mutant whose emergence we wish to prevent. If we are concerned about a specific VoC, then n = 1. If we are concerned about the set of all possible VoCs, then n>1, but n<L−m, where L is the length of the genome. Lastly, if we are concerned about any viable mutant, then n = L−m. In the course of our analysis, we found that the value of n has little effect on the evolutionary safety of a mutagenic treatment. As common in mutagenesis and also in the specific mechanism of action of Molnupiravir, transition mutations are more likely than transversion mutations (see Fig 1A). Our model can be extended to consider situations where the mutagenic drug increases the probability of mutation for a subset of all possible mutations (see Methods). Both x and y are cleared at same rate a j with the subscript j indicating the presence or absence of an adaptive immune response, such that during the growth phase j = 0 and during the clearance phase j = 1. We have a 0 <b<a 1 . Virus dynamics [16] in an infected patient can be described by the system of differential equations (1) We ignore back mutation from mutant to wild type [4,7,16]. In the growth phase, without treatment, we have bqm+n>a 0 since both x and y grow exponentially. In the clearance phase, without treatment, we have bqm<a 1 since both x and y decline exponentially. The system is linear and can be solved analytically (see Methods). The biological reactions are presented schematically in Fig 1B. In our simple approach, there is a sharp onset of adaptive immunity that happens at time T. We relax this assumption in a model extension (see section “Gradual activation of the immune system”).

Mutation rates We denote by u 0 the mutation rate without mutagenic treatment and by u 1 , which is greater than u 0 , the mutation rate with mutagenic treatment. The typical mutation rate for other positive single-strand RNA viruses is 10−5 per base [54]. The mutation rate of SARS‑CoV‑2 has been hypothesized to be lower because of a proofreading capability [55]. The per-base mutation rate has been estimated at u 0 = 10−6 per bp by proxy with the related beta-coronavirus MHV [56,57]. An in vitro study of experimental evolution of SARS‑CoV‑2 has reached the estimate u 0 = 3.7∙10−6 per base [58]. Zhou and colleagues [14] estimated the mutation rate in vitro of SARS‑CoV‑2 to be closer to 10−5 per base. Although the mode of replication (i.e., stamping versus linear replication [59]) and the number of mutations in the plus/minus strands can affect the distribution of the number of mutants, we found that it does not affect their expected number. For a detailed analysis of the influence of the mode of replication, the mutation rate in the plus/minus strand and RNA editing on the mutation rate, see section in File S3 “Estimating the mutation rate.” For our main analysis, we use u 0 = 10−6 per base. We also explore results for u 0 = 5∙10−6 and u 0 = 10−5. The mutation rate of SARS‑CoV‑2 under Molnupiravir treatment has been measured in vitro to be 2- to 5-fold higher than without treatment [14]. The fold increase in mutation rate under treatment can also be estimated from sequencing viral samples from treated patients. A 2-fold increase in the mutation rate in RNA-dependent RNA polymerase sequence in patients treated with Molnupiravir has been observed during its Phase 2a clinical trial [44]. This estimate comes with the caveat of neglecting potentially rare, severely deleterious mutants since those are less likely to be sequenced. Hence, we estimate u 1 to be 2 to 5 times higher than u 0 , relying on the in vitro estimate. Mutation rate estimations for different pathogen–drug combinations are available in the literature and result in even higher estimates for the virus mutation rate under treatment [20]. In our analysis, we explore a wide range of u 1 values, because it is our expectation that future mutagenic treatments might achieve higher increases of the virus mutation rate.

Viral birth and clearance rates The average lifetime of SARS‑CoV‑2 has been measured by proxy with MHV in monkey kidney cells and was found to be about 8 hours [57]. Hence, without infection of new cells, we would obtain a clearance rate of a 0 = 3 per day. From the current literature, we know that the virus load grows by about 10 orders of magnitude within 5 days [49,60]. Hence, for the viral growth rate we obtain b = 7.61 per virion per day. For the clearance phase, a decrease by 4 orders of magnitude in 10 days results in a death rate of a 1 = 8.76 per day reflecting high immunocompetence. The same fold decrease over 120 days results in a death rate of a 1 = 7.69 per day reflecting low immunocompetence (see Methods). These estimates are approximations as they ignore loss by lethal mutants.

Number of viral genome positions that are either lethal or potentially concerning when mutated The fitness effects of some individual mutations has been studied for some viruses such as the influenza virus [61,62], the HIV [63], or hepatitis C virus [64]. However, the number of lethal mutations m and the number of mutations that are potentially concerning needs to be computed for the whole SARS‑CoV‑2 genome. The distribution of fitness effects of random, single mutations has been studied in a different single-stranded RNA virus, the vesicular stomatitis virus (VSV) [65]. This distribution seems to be similar among single-stranded RNA viruses but could differ between species [66]. According to these studies, the proportion of mutations that are lethal when mutated is about 40% and the proportion of highly deleterious mutations, defined as those that reduce the viral fitness by more than 25%, represents about 30% of possible mutations. Note that the low mutation rate that is characteristic of SARS‑CoV‑2 allows us to approximate the number of lethal positions as 1/3 of the total number of possible mutations, taking into account that each position can be mutated to 3 different destinations (see Methods). SARS‑CoV‑2 genome has a length of 29,900 nucleotides. Hence, assuming 40% of positions being lethal upon mutation, we have m = 11,960 positions when considering lethal mutations only and m = 20,930 positions when considering that 70% of positions are either lethal or highly deleterious upon mutations. Hence, the realistic range for m is between 11,960 and 20,930 positions. For completeness, we also explore unrealistically lower bound of m such as 1,500 positions, which is the number of positions in the coding genome that are one nucleotide way from a STOP codon, assuming that most nonsense mutations are deleterious or lethal. We first consider mutants that are VoCs, which is exhibit phenotypes such as for example increased infectiousness or virulence. In order to estimate the number of positions that could give rise to new variants of concern when mutated (denoted by n), we used empirical data collected by [67,68]. Starr and colleagues conducted deep mutagenesis scans of the receptor-binding domain of the SARS‑CoV‑2 spike protein. For each of the generated mutants, Starr and colleagues measured the mutant’s binding affinity to ACE2 that is the receptor used by SARS‑CoV‑2 to enter the human cell. In a subsequent study, Starr and colleagues also measured each mutant’s affinity to antibodies in order to assess the ability of each mutant to escape the adaptive immune response and antibody treatments. Both escape from antibody and increased affinity to ACE2 are phenotypes beneficial for SARS‑CoV‑2 and are hence concerning. We identified 484 amino acid substitutions that result in antibody escape and 314 distinct amino acid substitutions that result in increased binding to ACE2. For each position coding for the receptor-binding domain of the spike protein, we counted how many mutations can give rise to the identified set of substitutions with either increased binding to ACE2 or decreased binding to antibodies (we corrected for the overlap of substitutions found in both categories). We found that the resulting estimate (divided by 3 to take into account all possible destinations, see Methods) was n = 87 positions when considering all possible mutations and n = 75 positions when considering only transition mutations, i.e., when taking into account the specific mechanism of action of Molnupiravir. Of course, mutations that are advantageous for the virus could occur also outside of the receptor-binding domain of the spike protein. More broadly, any neutral and even slightly deleterious mutation can be undesirable since they could represent an evolutionary “stepping-stone” to a multiple-mutation variant due to epistasis. Hence, we also explore how considering a very large number of positions that could be potentially concerning when mutated, up to the length of the SARS‑CoV‑2 genome minus the m positions that are lethal when mutated. Thus, our upper bound on the value of n is 29,900−m. In addition, we explore the possibility of double mutants in the section “Evolutionary safety for higher-order mutants.”

Abundance of mutant virus for various treatment regimes In Fig A1 in S1 Text, we show the dynamics of total virus and mutant over the course of an infection. We consider 4 times for the start of mutagenic treatment: at infection; at day 2 after infection, which corresponds to the beginning of symptoms; at day 5 after infection, which corresponds to the peak of the virus load; and at day 7 after infection. We observe that under each of these 4 options, treatment decreases the abundance of wild-type virus. The dynamics of mutant follows that of the wild type. For the parameters used in Fig A1 in S1 Text, treatment decreases the abundance of mutant virus—with exception of a brief transient period soon after the start of therapy, which is almost invisible in the figure. In Fig A2 in S1 Text, we assess the plausibility of the model by plotting virus load versus time for different values of the death rate during the clearance phase and comparing it to sequential measurements of virus load in patients. In Table A2 in S1 Text, we use measurements of virus load from patients that were treated or not with Molnupiravir. We are interested in calculating the total number of mutant virus produced over the course of infection. This number can be computed as the integral of the abundance of mutant virus over time (see Methods). We consider 2 scenarios: in the first, the patient begins treatment when their virus load reaches its peak; in the second, the patient begins treatment when they become infected (following exposure to an infected individual). Note that even without mutagenic treatment, due to the innate mutation rate, viral mutants will appear. Thus, our aim is to evaluate their total abundance for various mutation rates, with and without treatment.

Treatment begins at (or near) peak virus load In Fig 2, we show the cumulative mutant load, Y(u 1 ), as a function of the mutation rate u 1 for the case where treatment starts at peak virus load. To understand this function, we introduce the parameter η = x T /y T , with x T and y T denoting, respectively, wild-type and mutant virus load at peak. If η>n/m then Y(u 1 ) is a declining function. In this case, any mutagenic treatment is evolutionarily safe in the sense of reducing the cumulative mutant virus load compared to no treatment. If η<n/m then the function Y(u 1 ) attains a single maximum at (2) If u 1 >u* then any increase in mutation rate is beneficial for evolutionary safety as it actually decreases the chance of appearance of potentially concerning mutants compared to evolution of the virus under no treatment. If u 1 <u* then a small increase in the mutation rate can increase the chance of appearance of potentially concerning mutants under treatment, and thus be evolutionarily unsafe; in this case, there needs to be a sufficiently large increase in mutation rate to make the treatment evolutionarily safe (see Fig 2 for details). We notice that increasing estimates of m or decreasing a 1 reduces the value of u* and therefore increases the range of u 0 for which mutagenic treatment is evolutionarily safe. In particular, the slower the patient clears the virus (lower a 1 ), the lower the value of u* and hence treatments become more evolutionarily safe. In Fig 2, we notice that only for low m and high a 1 we find u*>u 0 . For all other cases, u*<u 0 , and lethal mutagenesis is both evolutionarily safe and desired, because it reduces the abundance of both wild type and mutant. PPT PowerPoint slide

PNG larger image

TIFF original image Download: Fig 2. Cumulative mutant virus load versus mutation rate, u 1 , during treatment. The cumulative mutant virus load increases with mutation rate u 1 before reaching a peak and then decreases to low values. If the peak is reached at a mutation rate that is less than the natural mutation rate, u 0 (red dotted line), then any increase in mutation rate reduces the cumulative mutant load. If the peak is reached for a mutation rate greater than u 0 , then the increase in mutation rate caused by mutagenic treatment must exceed a threshold value (blue dotted line) to reduce the cumulative mutant virus load. We also consider mutants with a 1% advantage in the birth rate. As expected, we observe a higher cumulative mutant load for the advantageous mutant (green line) compared to the neutral mutant (blue line). But the minimum mutation rate under treatment that is required for evolutionary safety is slightly lower for the advantageous mutant. (A) Treatment starts at peak virus load. (B) Treatment starts at infection. The red arrow indicates the mutation rate at the error threshold of the growth phase. Parameters: b = 7.61 per day, a 0 = 3 per day, n = 1 position, T = 5 days, m and a 1 as shown. The code used to generate this figure can be found at DOI: 10.5281/zenodo.8017992. https://doi.org/10.1371/journal.pbio.3002214.g002 If the mutagenic treatment is strong enough, it precludes the replication of the virus. Hence, it is always safe. We indicate, with a red arrow, the error threshold for the virus (see Fig 2). If the mutation rate exceeds the error threshold, the probability of lethal mutation occurring with each viral replication is so high that the viral population cannot grow, and hence, the population becomes extinct [16]. The error threshold is the mutation rate u e that solves b (1−u e )m = a 0 .

Treatment begins at (or soon after) infection In Fig 2, we also show the cumulative mutant load, Y(u 1 ), as a function of the mutation rate u 1 for the case where treatment starts at infection. We find that this function attains a maximum at a value which is given by the root of a third order polynomial (see Methods and Fig A3 in S1 Text). Using the notation and h = bT, we can approximate u* as follows: (3) Again if u 0 >u* then any increase in mutation rate is beneficial. If u 0 <u* then a small increase in the mutation rate can be evolutionarily not safe, but a sufficiently large increase in mutation rate can make the treatment evolutionarily safe (see Fig 2 for more details). We also find that early treatment is more effective in reducing mutant virus load when compared to late treatment.

Exploring the parameter space for evolutionary safety In Fig 3, we show the fold increase in virus mutation rate that mutagenic treatment has to achieve beyond the innate mutation rate of the virus to be evolutionarily safe. We vary first the estimated number of lethal mutations m in the viral genome and the clearance rate a 1 . For treatment starting at peak virus load (Fig 3A), we find that increase in mutation rate is evolutionarily safe if m>22,000 or a 1 <7.8 per day (green region). Evolutionary safety becomes an issue for small values of m and larger values of a 1 . For m = 12,000 positions and a 1 = 9 per day, we need at least a 10-fold increase in mutation rate before the drug attains evolutionary safety. When treatment begins at infection (Fig 3B), the evolutionarily safe area becomes smaller, but the minimum increase in mutation rate required for evolutionary safety is lower. For example, for a 1 = 9 per day and m = 12,000 positions, we need only a 3-fold increase. We show the same figure, but for an extended range of m values in Fig A4 in S1 Text. Note that our estimate for fold increase in mutation rate for Molnupiravir is about 2, so the drug is safe only for a portion of the parameter space. PPT PowerPoint slide

PNG larger image

TIFF original image Download: Fig 3. Evolutionary safety of mutagenic treatment. In the green parameter region, any increase in mutation rate reduces the cumulative mutant virus load and is therefore evolutionarily safe. In the red shaded region, we indicate the minimum fold increase in mutation rate that is required to reduce the cumulative mutant load. Contour lines for 3-fold and 10-fold increase are shown. (A) Treatment starts at peak virus load. (B) Treatment starts at infection. Parameters: b = 7.61 per day, a 0 = 3 per day, n = 1, T = 5 days, u 0 = 10−6 per bp. The code used to generate this figure can be found at DOI: 10.5281/zenodo.8017992. https://doi.org/10.1371/journal.pbio.3002214.g003

The evolutionary risk factor is a slowly declining function of the number of mutations leading to viable virus So far, we have used the parameter n = 87 to denote the number of mutations that would result in VoCs that is variants with increased transmissibility, virulence, or resistance to existing vaccines and treatments. However, in the broad sense, any treatment that increases the standing genetic variation of the virus could favor the emergence of new variants of concern by enabling epistatic mutations. Therefore, we now extend the interpretation of n to include any viable mutation in the viral genome. In Fig 5, we show that the ERF is a declining function of n. Thus, the more opportunities the virus has for viable mutations (the larger n), the higher the advantage of mutagenic treatment. The reason for this counter-intuitive observation is that for large n the cumulative mutant virus load is high already in the absence of treatment, while mutagenic treatment reduces the mutant load by forcing additional lethal mutations. ERF decreases with the number of positions n also for lower birth rate b (Fig A13 in S1 Text). PPT PowerPoint slide

PNG larger image

TIFF original image Download: Fig 5. The ERF versus the number, n, of positions in the viral genome giving rise to concerning (or viable) mutations. The ERF of mutagenic treatment is the ratio of the cumulative mutant virus load with and without treatment. We explore all values of n subject to the constraint that m+n remains below the length of the SARS‑CoV‑2 genome. We observe that the ERF decreases as function of n. (A) Treatment starts at peak virus load. (B) Treatment starts at infection. Parameters: a 0 = 3 per day, b = 7.61 per day, u 0 = 10−6 per bp, u 1 = 3∙10−6 per bp, T = 5 days. Initial condition: x 0 = 1 and y 0 = 0. The code used to generate this figure can be found at DOI: 10.5281/zenodo.8017992. ERF, evolutionary risk factor; SARS‑CoV‑2, Severe Acute Respiratory Syndrome Coronavirus 2. https://doi.org/10.1371/journal.pbio.3002214.g005

Advantageous mutants do not substantially affect the evolutionary safety compared to neutral mutants Mutants could have an in-host advantage compared to wild type, such as faster a reproductive rate or a lower clearance rate. In Fig 2, we evaluate a mutant with a 1% selective advantage in birth rate. We have also included results for a mutant with a 0.5% selective advantage in Fig A14 in S1 Text, and results for mutants with 0.5% and 1% selective disadvantage in Figs A15 and A16 in S1 Text. As expected, we observe that the advantageous mutant reaches higher virus load than a neutral mutant. But we also observe that if there is a minimum increase in mutation rate that is required for evolutionary safety, then it is lower (or slightly lower) for the advantageous mutant. Therefore, a treatment that is evolutionarily safe for a neutral mutant is also evolutionarily safe for an advantageous mutant. Conversely, we observe that the minimum increase in the mutation rate required for evolutionary safety is higher when the mutant has a selective disadvantage. We also consider a 1% disadvantage of the wild type with regards to the mutant that is exhibited under treatment, that is when u 0 >10−6. In this scenario, we find that the treatment is always evolutionarily safe (see Fig A17 in S1 Text).

Gradual activation of the immune system So far, we have considered a sudden activation of the adaptive immune response by switching the clearance from a 0 to a 1 at time T resulting in a two-phase model of immunity. In reality, the immune response intensifies gradually over the course of the infection [16]. We explore a more gradual onset of the immune response in Fig A18 in S1 Text, where we add an intermediate phase during which the clearance rate is the arithmetic average of a 0 and a 1 . We find that the ERF value for the three-phase immunity is very close to and bounded by the ERF values found for corresponding two-phase simulations.

Nonlethal deleterious mutations In order to make the model more biologically realistic, we have extended our model to consider nonlethal deleterious mutations. In our extended model, we consider 4 categories of virus: the wild-type; nonlethal deleterious mutants with no concerning mutations; mutants with concerning mutations; and nonlethal deleterious mutant with concerning mutations. We found that considering nonlethal deleterious mutations always increases the evolutionary safety of the treatment. For a detailed analysis, see section “Non-lethal deleterious mutants” in S3 Text.

Lethal defection Nonviable virions can negatively interfere with the replication of viable virions that are being generated from the same cell. For example, nonviable virions may consume functional proteins that are then lacking for the replication and packaging of the viable virions. Conversely, nonfunctional proteins synthetized from the mutant genomes can lead to the production of noninfective virus particles containing wild-type genomes. We extended our model to take into account the interference of nonviable mutants in the replication of viable mutants. We found that the evolutionary safety is always increased when incorporating this effect. For a detailed analysis, see S3 Text.

A simple approach captures the essence of mutagenic treatment and evolutionary safety We further simplify our mathematical framework to obtain quantitative guidelines about the evolutionary safety of a mutagenic drug. We find that focusing on virus dynamics in the growth phase can be used to approximate the full infection dynamics, especially if the clearance rate is large. Note that clearance rates leading to infections which last longer than 100 days remain exceptions, and hence, most individuals have a high clearance rate a 1 . The simplified approach is presented in the Methods. The agreement between the simplified and the full model is shown in Fig A19 in S1 Text. The eventual goal of all mutagenic treatments would be to prevent the exponential expansion of the virus even before the onset of adaptive immunity. Using the SARS‑CoV‑2 estimates, m = 20,000 positions, b = 7.61 per day, and a 0 = 3 per day, we find that mutagenic treatment would have to achieve u 1 >4.65∙10−5, which is a 50-fold increase of the natural mutation rate of the virus. If the mutagenic drug results in a smaller increase in the virus mutation rate under treatment, then it does not prevent the establishment of the infection, but it could still reduce both wild type and mutant abundance. The mutant virus load at time T is a one-humped function of the mutation rate with a maximum that is close to u* = 1/(bTm). For m = 20,000 positions, b = 7.61 per day, and T = 5 days, we find u* = 1.32∙10−6 per bp. This value is close to the estimate of the natural mutation rate of the virus, u 0 = 10−6 per bp. If u 0 was greater than u* then any increase in mutation rate would be evolutionarily safe. Otherwise, we need to calculate the condition for evolutionary safety. Let us introduce the parameter s with u 1 = su 0 . The condition for evolutionary safety in the simplified model is (4) As before b = 7.61 per day, T = 5 days, and u 0 = 10−6 per bp. For s = 3 fold increase of mutation rate induced by mutagenic treatment, we get m>14,455 positions. Since evolutionary safety improves with decreasing clearance rate a 1 (in the full model), we can interpret inequality (4) as a sufficient condition or as an upper bound. The agreement between the analytical formulas and the numerical computation of the model is shown in Fig A20 in S1 Text. For the simplified model, we also find that ERF is a declining function of the number of viable mutations, n (see Fig A21 in S1 Text).

Evolutionary safety for higher-order mutants Up until now, we considered only one-step mutations in order to produce viable virus. However, many examples of multistep adaptation have been observed in virus evolution. In particular, the H274Y mutation has been studied and found to be deleterious, but the loss in fitness could be restored by other mutations [61,62]. In this section, we consider double mutants. The evolutionary dynamics for a virus with a two-locus, binary genome can be written as: (5) The frequencies y 00 , y 01 , y 10 , and y 11 denote, respectively, the wild-type virus (00), 2 single mutants (01 and 10), and the double mutant (11). We have q = 1−u, where q is the quality of replication, and u is the mutation rate. We have B ij = (1+ds ij )b, where B ij the birth rate of the variant ij, d is the fitness difference between the wild type and the mutant, and s ij determines whether the mutant has a fitness disadvantage (s ij = −1), is neutral (s ij = 0) or has a fitness advantage (s ij = 1). For the wild type (00), we always set s 00 = 0, hence B 00 = b. The fitness landscape is given by the vector (s 01 , s 10 , s 11 ). Since each component has one of 3 values (−1, 0, or 1), there are 27 possible landscapes. Because of symmetry, 9 landscapes are redundant. Hence, we are left with 18 different landscapes. When considering double mutants, we propose 2 different methods for evaluating evolutionary safety. In Method 1, evolutionary safety requires that treatment reduces the sum of all mutants (single and double). In this case, double mutants make only a negligible contribution to evolutionary safety because their abundance is much lower than that of single mutants (see Figs A22–A24 in S1 Text). In Method 2, evolutionary safety requires that treatment reduces the amount of each mutant type separately. Here, we find that in some cases evolutionary safety requires a larger increase in mutation rate in order to reduce the amount of double mutant (see Figs A25–A27 in S1 Text). In Figs A28 and A29 in S1 Text as well as in Tables A3 and A4 in S1 Text, we investigate the effect of different fitness landscapes on evolutionary safety. In Fig A26 in S1 Text and in Table A3 in S1 Text, we consider a patient that has a low clearance rate of infection. In Fig A29 in S1 Text and in Table A4 in S1 Text, we consider a patient that has a high clearance rate of infection (where in general mutagenic treatment is less likely to be evolutionarily safe). In Fig A28 in S1 Text, we see that any increase in mutation rate reduces the cumulative sum of all mutants and is therefore evolutionarily safe using Method 1. For 4 of the 5 landscapes shown here, any increase in mutation rate reduces the abundance of the double mutant and is therefore evolutionary safe using Method 2. In Fig A29 in S1 Text, we see that any increase in mutation rate reduces the cumulative sum of all mutants and is therefore evolutionarily safe using Method 1. For all 5 landscapes shown here, we need roughly a 3-fold increase in mutation rate to reduce the abundance of the double mutant and therefore achieve evolutionary safe using Method 2. In Table A3 in S1 Text, we see that for 14 of the 18 fitness landscapes, any increase in mutation rate reduces the abundance of the double mutant, but for 4 landscapes we need an increase in mutation rate between 1.2- and 1.5-fold to reduce the abundance of the double mutant. In Table A4 in S1 Text, we see that for all 18 fitness landscapes, we need approximately a 3.5- to 3.8-fold increase in mutation rate to reduce the abundance of the double mutant.

[END]
---
[1] Url: https://journals.plos.org/plosbiology/article?id=10.1371/journal.pbio.3002214

Published and (C) by PLOS One
Content appears here under this condition or license: Creative Commons - Attribution BY 4.0.

via Magical.Fish Gopher News Feeds:
gopher://magical.fish/1/feeds/news/plosone/