Abstract

We investigated the mutation profiles of severe acute respiratory syndrome coronavirus 2 in samples collected from a molnupiravir and nirmatrelvir/ritonavir combination therapy in macaques. We found that molnupiravir induced several nirmatrelvir resistance mutations at low abundance that were not further selected in combination therapy. Coadministration of nirmatrelvir/ritonavir lowered the magnitude of the mutagenetic effect of molnupiravir.

The severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) pandemic has persisted for >4 years, presenting an ongoing risk to public health. Current treatments for early symptomatic SARS-CoV-2 rely on several direct-acting antivirals (DAAs). Three small-molecule compounds have been authorized for treating SARS-CoV-2 infection in high-risk populations: the 3C-like protease (3CLpro) inhibitor nirmatrelvir/ritonavir (NMV/r), the mutagenic nucleoside analog molnupiravir (MOV), and the RNA-dependent RNA polymerase (RdRp) chain-terminator remdesivir [1]. There have been concerns raised regarding MOV that relate to its mutagenic mechanism. The most obvious concern is the potential for incorporation of MOV into the host DNA, leading to genetic mutations in host cells [2, 3]. Additionally, there are concerns that increased viral genome mutations could enhance the likelihood of the generating new variants [4]. These variants could include resistance mutations against other DAAs such as NMV/r, or antibody escape mutations that would have reduced control by immunity from earlier infection and vaccination.

The present treatment approach for SARS-CoV-2 primarily relies on a 5-day course of a single DAA, with NMV/r being the most frequently prescribed DAA in the United States [5]. The combination of the oral DAAs NMV/r and MOV may provide addition therapeutic benefit, but this approach has not been tested in clinical trials. In a previous study, we compared the antiviral effect of MOV or NMV/r and the coadministration of both in a SARS-CoV-2 macaque model. Macaques received the oral treatment of MOV, NMV/r, combination therapy, or vehicle control 12 hours after exposure to SARS-CoV-2, and then were followed for 4 days before necropsy. The combination therapy resulted in milder disease progression, a greater reduction in virus titer, and reduced lung pathology compared to either single-agent treatment [6].

The earlier study offered a unique opportunity to study the combined effects of MOV and NMV/r and to examine the influence of the mutagenic effects of MOV on the emergence of resistance mutations to NMV/r. Specifically, we used a highly accurate viral RNA sequencing approach, multiplexed Primer ID next-generation sequencing (MPID-NGS) [7], to examine the SARS-CoV-2 mutation profiles from the different treatment groups. We found that NMV/r reduced MOV-induced mutagenicity of SARS-CoV-2 when coadministered, while still contributing benefit to dual therapy [6]. The mutagenic property of MOV increased several mutations related to 3CLpro resistance, but these mutations were not further enriched in the combination therapy.

METHODS

We used a previously published MPID-NGS protocol to deep-sequence multiple regions of the viral genome of SARS-CoV-2 [2]. Total tissue RNA was extracted from tissue samples taken at necropsy from each lobe of the lungs in the macaques after 4 days of treatment [6]. The lung lobe with the highest SARS-CoV-2 titer in each animal, measured using quantitative polymerase chain reaction, was used for the library preparation. The sequenced regions included all of nsp5 (3CLpro), a portion of nsp12 (RdRp), and the S gene receptor-binding domain (S-RBD). Information regarding the approval for the animal study can be found in the Supplementary Materials. Sequences of the primers used in the library prep are described in Supplementary Table 1.

The sequence data were processed using the tcs pipeline (v.2.5.2) to construct a template consensus sequence for each viral template/genome sequenced. We then calculated the mutation density at each sequenced region. We distinguished true substitutions from the residual background method error (0.01%) by calculating a false discovery rate (FDR)–adjusted P value [8]. Substitutions with an FDR-adjusted P < .05 were considered as true mutations in the viral population. The statistical analysis and graphic illustration were performed using R software (version 4.1.2).

RESULTS

Mutation Densities in Combination Therapy Versus Single DAA Therapy

We sequenced SARS-CoV-2 viral RNA extracted from lung tissues obtained from a total of 20 animals, divided into 4 treatment groups (MOV, NMV/r, MOV + NMV/r, or vehicle only), with each group containing 5 animals. Sequencing was not successful in 2 animals from the MOV + NMV/r group because of low viral RNA titers in the tissue, consistent with the additive effect of these inhibitors on control of virus replication and clinical outcome [6]. The median number of template consensus sequences for each animal, which represents the number of viral RNA templates/genomes sequenced at the 3 regions, were as follows: 2888 at nsp5, 2207 at nsp12, and 3004 at the S-RBD region.

We first examined the mutation density in SARS-CoV-2 RNA from the 4 groups of animals at the 3 sequenced regions (Figure 1A). In a previous study, we showed that the active form of MOV, β-d-N4-hydroxycytidine (NHC), preferentially pairs with G (as C), and is sometimes read as a U by polymerase, causing C-to-U or G-to-A mutations [2, 9]. In this study, our results were consistent with this previous finding. In the MOV group, the overall substitution rate was at approximately 0.05%, significantly higher than other groups, driven by the increase of the C-to-U and G-to-A mutations. The overall substitution rates were significantly increased in the MOV + NMV/r group compared with NMV/r or the vehicle groups, but lower than the MOV group. The C-to-U mutation and G-to-A mutation densities in the combination therapy group followed the same trend as seen in the MOV monotherapy. These results suggest that the coadministration of NMV/r lowered the mutagenetic effect of MOV against SARS-CoV-2.

Measurement of mutation rates and assessment of potential resistance mutations on 3C-like protease (3CLpro). A, Mutation density of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) from 4 groups of animals (molnupiravir [MOV], nirmatrelvir/ritonavir [NMV/r], MOV + NMV/r, and vehicle, shown in this order left to right) at 3 sequenced regions. Sequencing was performed using multiplexed Primer ID next-generation sequencing at nsp5 (3CLpro) and nsp12 (RNA-dependent RNA polymerase) and S gene receptor-binding domain regions of SARS-CoV-2. Substitution rates between groups were compared using t test. Significance levels: ***P < .001; **P < .01; *P < .05. B, Mutation density at each nsp5/3CLpro amino acid positions. C, Potential 3CLpro drug resistance mutations against NMV/r, mutation rate per 10 000 amino acids. Mutations with false discovery rate–adjusted P < .05 are highlighted. Abbreviations: 3CLpro, 3C-like protease; AA, amino acid; EC50, half maximal effective concentration; MOV, molnupiravir; ND, no data; NMV/r, nirmatrelvir/ritonavir; nts, nucleotides; S-RBD, S gene receptor-binding domain.
Figure 1.

Measurement of mutation rates and assessment of potential resistance mutations on 3C-like protease (3CLpro). A, Mutation density of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) from 4 groups of animals (molnupiravir [MOV], nirmatrelvir/ritonavir [NMV/r], MOV + NMV/r, and vehicle, shown in this order left to right) at 3 sequenced regions. Sequencing was performed using multiplexed Primer ID next-generation sequencing at nsp5 (3CLpro) and nsp12 (RNA-dependent RNA polymerase) and S gene receptor-binding domain regions of SARS-CoV-2. Substitution rates between groups were compared using t test. Significance levels: ***P < .001; **P < .01; *P < .05. B, Mutation density at each nsp5/3CLpro amino acid positions. C, Potential 3CLpro drug resistance mutations against NMV/r, mutation rate per 10 000 amino acids. Mutations with false discovery rate–adjusted P < .05 are highlighted. Abbreviations: 3CLpro, 3C-like protease; AA, amino acid; EC50, half maximal effective concentration; MOV, molnupiravir; ND, no data; NMV/r, nirmatrelvir/ritonavir; nts, nucleotides; S-RBD, S gene receptor-binding domain.

Detection of Reported Nirmatrelvir Resistance Mutations

We performed 2 analyses to study the impact of MOV on the potential development of resistance mutations against NMV/r in nsp5/3CLpro. We first plotted the substitution rate at each amino acid position of nsp5/3CLpro in each group (Figure 1B). The goal of this analysis was to visualize any potential mutation hotspots on the region that could be amplified by selection as a resistance mutation. The overall mutation rate was higher in the MOV group and the combination group, but we did not observe any mutation hotspots. In the second analysis, we assessed the amino acid codon changes at 8 positions (T21I, L50F, E166A, E166V, L167F, A173V, P252L, and T304I) that have been reported to contribute to NMV resistance, along with their reported half maximal effective concentration value fold-change in in vitro cell culture (Figure 1C) [10]. Mutations with FDR-adjusted P < .05, representing true mutations above background, are highlighted. Several resistance mutations including T21I, L50F, A173V, and P252L were significantly increased in the MOV group but were not further enriched in the combination treatment group. All of these mutations were associated with C-to-U mutations. Although the overall intrahost frequencies of these mutations were low (∼0.1%–0.2%), they were significantly higher than the background error level, but not higher than the average C-to-U MOV-induced mutation rate across the nsp5 region (Figure 1B). Several other NMV/r-related resistance mutations, particularly E166V, were not significantly increased in any of the treatment groups. In summary, our findings suggest that MOV treatment alone increased the frequencies of several mutations related to NMV resistance involving C-to-U mutations, but there was no evidence that combination therapy increased the specific accumulation of these NMV resistance mutations through further selection.

DISCUSSION

In this study, we used specimens from a previous MOV and NMR/r combination SARS-CoV-2 treatment study in nonhuman primates to assess the potential effect of the 2 compounds on the accumulation of NMR/r resistance. We were specifically interested in the possibility of inducing NMR/r resistance mutations through mutagenicity of the coadministered MOV. To achieve this aim, we used an ultra-deep MPID-NGS viral sequencing approach that has an error rate as low as 0.01% [7].

Several studies have reported the mutation densities of coronaviruses with the treatment of NHC/MOV in in vitro and mice models, as well as in specimens collected from clinical trials [2, 9, 11, 12]. The mutation rate and mutation spectrum in the MOV group in this study agreed with these previous reports. The overall 0.05% mutation density induced by NHC/MOV represents approximately 15 random substitutions across the 30 kb SARS-CoV-2 genome. This number is high enough to ensure that essentially all genomes contain some number of random substitutions. Particularly, since the predominant substitution is C to T (C to U on the viral genome), this turns every CAA, CAG, and CGA mutation into a target for introducing a termination codon.

MOV relies on viral replication to induce mutations in SARS-CoV-2. However, in combination therapy, viral replication is directly inhibited by NMV/r due to its effect on 3CLpro and polyprotein cleavage, resulting in a reduced mutation rate in the viral genome. Even given the reduced mutation load induced by MOV in combination therapy, the fitness of the virus appears to be diminished, thereby enhancing the therapeutic effect beyond that achieved by MOV or NMV/r alone [6].

The absence of any increased frequency of specific mutations associated with NMV/r resistance in the combination therapy suggests it is unlikely that MOV treatment is potentiating the evolution of drug resistance against NMV/r when coadministered. Even though several mutations were significantly increased in the MOV group, variants containing the mutations were rarely seen in the viral population (∼0.1%–0.2%), and they were not higher than other positions that could be affected by C-to-U mutations. However, there is still a risk that these variants could be induced by suboptimal MOV treatment, potentially impacting future NMV/r treatment. A recent study reported that the introduction of MOV as a DAA against coronavirus disease 2019 was associated with increased MOV signature mutations found in countries with high usage of MOV [4]. Of note, the sequencing analysis of that study focused only on the MOV-associated G-to-A mutations, but in our animal study, we found that all of the NMV resistance mutations that were associated with MOV usage were C-to-U mutations.

In conclusion, our work supports the use of combination therapy as it had a greater therapeutic effect and did not increase the potential for NMV mutations. However, the MOV monotherapy–associated C-to-U mutations in viruses circulating in the human population and their potential impact on NMV resistance need further investigation.

Supplementary Data

Supplementary materials are available at The Journal of Infectious Diseases online (http://jid.oxfordjournals.org/). Supplementary materials consist of data provided by the author that are published to benefit the reader. The posted materials are not copyedited. The contents of all supplementary data are the sole responsibility of the authors. Questions or messages regarding errors should be addressed to the author.

Notes

Acknowledgments. We acknowledge the contribution of the University of North Carolina High Throughput Sequencing Facility.

Financial support. This work was supported by the National Institutes of Health Antiviral Drug Discovery and Development Center (grant number U19 AI142759); the National Institutes of Health (grant number R01 AI140970 to R. S.); and the Intramural Research Program of the National Institute of Allergy and Infectious Diseases. This research received infrastructure support from the University of North Carolina (UNC) Center for AIDS Research (grant number P30-AI050410) and the UNC Lineberger Comprehensive Cancer Center (grant number P30-CA016086).

References

1

Meyerowitz
 
EA
,
Li
 
Y
.
Review: the landscape of antiviral therapy for COVID-19 in the era of widespread population immunity and Omicron-lineage viruses
.
Clin Infect Dis
 
2023
;
78
:
908
17
.

2

Zhou
 
S
,
Hill
 
CS
,
Sarkar
 
S
, et al.   
β-d-N4-hydroxycytidine inhibits SARS-CoV-2 through lethal mutagenesis but is also mutagenic to mammalian cells
.
J Infect Dis
 
2021
;
224
:
415
9
.

3

Swanstrom
 
R
,
Schinazi
 
RF
.
Lethal mutagenesis as an antiviral strategy
.
Science
 
2022
;
375
:
497
8
.

4

Sanderson
 
T
,
Hisner
 
R
,
Donovan-Banfield
 
I
, et al.   
A molnupiravir-associated mutational signature in global SARS-CoV-2 genomes
.
Nature
 
2023
;
623
:
594
600
.

5

Murphy
 
SJ
,
Samson
 
LW
,
Sommers
 
BD
.
COVID-19 antivirals utilization: geographic and demographic patterns of treatment in 2022
.
Washington, D.C.
:
Office of the Assistant Secretary for Planning and Evaluation, US Department of Health and Human Services
,
2022
.

6

Rosenke
 
K
,
Lewis
 
MC
,
Feldmann
 
F
, et al.   
Combined molnupiravir-nirmatrelvir treatment improves the inhibitory effect on SARS-CoV-2 in macaques
.
JCI Insight
 
2023
;
8
:
e166485
.

7

Zhou
 
S
,
Jones
 
C
,
Mieczkowski
 
P
,
Swanstrom
 
R
.
Primer ID validates template sampling depth and greatly reduces the error rate of next-generation sequencing of HIV-1 genomic RNA populations
.
J Virol
 
2015
;
89
:
8540
55
.

8

Zhou
 
S
,
Hill
 
CS
,
Spielvogel
 
E
,
Clark
 
MU
,
Hudgens
 
MG
,
Swanstrom
 
R
.
Unique molecular identifiers and multiplexing amplicons maximize the utility of deep sequencing to critically assess population diversity in RNA viruses
.
ACS Infect Dis
 
2022
;
8
:
2505
14
.

9

Sheahan
 
TP
,
Sims
 
AC
,
Zhou
 
S
, et al.   
An orally bioavailable broad-spectrum antiviral inhibits SARS-CoV-2 in human airway epithelial cell cultures and multiple coronaviruses in mice
.
Sci Transl Med
 
2020
;
12
:
eabb5883
.

10

US Food and Drug Administration. Fact sheet for healthcare providers: emergency use authorization for Paxlovid. 2021. https://www.fda.gov/media/155050/download. Accessed 20 November 2023.

11

Agostini
 
ML
,
Pruijssers
 
AJ
,
Chappell
 
JD
, et al.   
Small-molecule antiviral β-d-N(4)-hydroxycytidine inhibits a proofreading-intact coronavirus with a high genetic barrier to resistance
.
J Virol
 
2019
;
93
:
e01348-19
.

12

Fischer
 
WA
,
Eron
 
JJ
,
Holman
 
W
, et al.   
A phase 2a clinical trial of molnupiravir in patients with COVID-19 shows accelerated SARS-CoV-2 RNA clearance and elimination of infectious virus
.
Sci Transl Med
 
2022
;
14
:
eabl7430
.

Author notes

Potential conflicts of interest. The University of North Carolina is pursuing intellectual property protection for Primer ID and R. S. has received nominal royalties. All other authors report no potential conflicts.

All authors have submitted the ICMJE Form for Disclosure of Potential Conflicts of Interest. Conflicts that the editors consider relevant to the content of the manuscript have been disclosed.

This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/pages/standard-publication-reuse-rights)

Supplementary data