Genomics and epidemiology of a novel SARS-CoV-2 lineage in Manaus, Brazil.

Cases of SARS-CoV-2 infection in Manaus, Brazil, resurged in late 2020, despite high levels of previous infection there. Through genome sequencing of viruses sampled in Manaus between November 2020 and January 2021, we identified the emergence and circulation of a novel SARS-CoV-2 variant of concern, lineage P.1, that acquired 17 mutations, including a trio in the spike protein (K417T, E484K and N501Y) associated with increased binding to the human ACE2 receptor. Molecular clock analysis shows that P.1 emergence occurred around early November 2020 and was preceded by a period of faster molecular evolution. Using a two-category dynamical model that integrates genomic and mortality data, we estimate that P.1 may be 1.4-2.2 times more transmissible and 25-61% more likely to evade protective immunity elicited by previous infection with non-P.1 lineages. Enhanced global genomic surveillance of variants of concern, which may exhibit increased transmissibility and/or immune evasion, is critical to accelerate pandemic responsiveness. ONE-SENTENCE
SUMMARY: We report the evolution and emergence of a SARS-CoV-2 lineage of concern associated with rapid transmission in Manaus.

Brazil has experienced high mortality during the COVID-19 pandemic, recording >250,000 deaths and >10 million reported cases, as of February 2020. SARS-CoV-2 infection and disease burden have been highly variable across the country, with Amazonas state in north Brazil being the worst-affected region (1). Serological surveillance of blood donors in Manaus, the capital city of Amazonas and the largest city in the Amazon region, has suggested >67% cumulative attack rates by October 2020 (2). Similar but slightly lower seroprevalences have also been reported for cities in neighbouring regions (3, 4). However, the level of previous infection in Manaus was clearly not sufficient to prevent a rapid resurgence in SARS-CoV-2 transmission and mortality there during late 2020 and early 2021 (5), which has placed a significant pressure on the city’s healthcare system.

Here, we show that the second wave of infection in Manaus was associated with the emergence and rapid spread of a new SARS-CoV-2 lineage of concern, named lineage P.1. The lineage carries a unique constellation of mutations, including several that have been previously determined to be of virological importance (6–10) and which are located in the spike protein receptor binding domain (RBD), the region of the virus involved in recognition of the angiotensin-converting enzyme-2 receptor cell surface receptor (11). Using genomic data, structure-based mapping of mutations of interest onto the spike protein, and dynamical epidemiology modelling of genomic and mortality data, we investigate the emergence of the P.1 lineage and explore epidemiological explanations for the resurgence of COVID-19 in Manaus.

Identification and nomenclature of a novel P.1 lineage in Manaus

In late 2020, two SARS-CoV-2 lineages of concern were discovered through genomic surveillance, both characterised by sets of significant mutations: lineage B.1.351, first reported in South Africa (12) and lineage B.1.1.7, detected in the United Kingdom (13). Both variants have transmitted rapidly in the countries where they were discovered and spread to other regions (14, 15). Analyses indicate B.1.1.7 has higher transmissibility than previously circulating lineages in the UK (13).

Following a rapid increase in hospitalizations in Manaus caused by severe acute respiratory infection in Dec 2020 (Fig. 1A), we focused ongoing SARS-CoV-2 genomic surveillance (2, 16–20) on recently collected samples from the city (see Materials and Methods). Prior to this, only seven SARS-CoV-2 genome sequences from Amazonas state were publicly available (SARS-CoV-2 was first detected in Manaus on 13 March 2020) (17, 21). We sequenced SARS-CoV-2 genomes from 184 samples from patients seeking COVID-19 testing in two diagnostic laboratories in Manaus between November and December 2020, using the ARTIC V3 multiplexed amplicon scheme (22) and the MinION sequencing platform. As partial genome sequences can provide useful epidemiological information, particularly regarding virus genetic diversity and lineage composition (23), we harnessed information from partial (n=41, 25–75% genome coverage), as well as near-complete (n=95, 75–95%) and complete (n=48, ≥95%) sequences from Manaus (fig. S1–S3), together with other available and published genomes from Brazil for context (Data S3). Pango lineages were classified using the Pangolin (24) software tool (http://pangolin.cog-uk.io/) and standard phylogenetic analysis using complete reference genomes.

Fig. 1.

SARS-CoV-2 epidemiological, diagnostic, genomic and mobility data from Manaus.

(A) Dark solid line shows the 7-day rolling average of the COVID-19 confirmed and suspected daily time series of hospitalisations in Manaus. Admissions in Manaus are from Fundação de Vigilância em Saúde do Amazonas (75). Green dots represent daily severe acute respiratory mortality records from the SIVEP-Gripe (Sistema de Informação de Vigilância Epidemiológica da Gripe) database (71). SARI = severe acute respiratory infections. Excess burial records based on data from Manaus Mayor’s office are shown in red dots for comparison (see Materials and Methods). The arrow denotes 6 December 2020, the date of the first P.1 case identified in Manaus by our study. (B) Maximum likelihood tree (n=974) with B.1.1.28, P.1 and P.2 sequences, with collapsed views of P.1 and P.2 clusters and highlighting other sequences from Manaus, Brazil). Ancestral branches leading to P.1 and P.2 are shown as dashed lines. See fig. S3 for a more detail phylogeny. Scale bar is shown in units of nucleotide substitutions per site (s/s). (C) Number of air travel passengers from Manaus to all states in Brazil was obtained from National Civil Aviation Agency of Brazil (www.gov.br/anac). The ISO 3166–2:BR codes of the states with genomic reports of P.1 (GISAID (76), as of 24 Feb 2021), are shown in bold. An updated list of GISAID genomes and reports of P.1 worldwide is available at https://covlineages.org/global_report_P.1.html. (D) Number of genome sequences from Manaus belonging to lineages of interest (see Materials and Methods); spike mutations of interest are denoted.

Our early data indicated the presence of a novel SARS-CoV-2 lineage in Manaus containing 17 amino acid changes (including 10 in the spike protein), 3 deletions, 4 synonymous mutations and a 4 nucleotide insertion compared to the most closely related available sequence (GISAID ID: EPI_ISL_722052) (25) (Fig. 1B). This lineage was given a new designation, P.1, on the basis that (i) it is phylogenetically and genetically distinct from ancestral viruses, (ii) associated with rapid spread in a new area, and (iii) carries a constellation of mutations that may have phenotypic relevance (24). Phylogenetic analysis indicated that P.1, and another lineage, P.2 (17), were descendants of lineage B.1.1.28 that was first detected in Brazil in early March 2020 (Fig 1B). Our preliminary results were shared with local teams on 10 Jan 2021 and published online on 12 Jan 2021 (25). Concurrently, cases of SARS-CoV-2 P.1 infection were reported in Japan in travellers from Amazonas (26). As of 24 Feb 2021, P.1 had been confirmed in 6 Brazilian states, which in total received >92,000 air passengers from Manaus in November 2020 (Fig. 1C). Genomic surveillance first detected lineage P.1 on 6 December 2020 (Fig. 1A), after which the frequency of P.1 relative to other lineage increased rapidly in the tested samples from Manaus (Fig. 1D; lineage frequency information can be found in fig. S2).

Dating the emergence of the P.1 lineage

We next sought to understand the emergence and evolution of lineage P.1 using molecular clock phylogenetics (23). We first regressed root-to-tip genetic distances against sequence sampling dates (27) for the P.1, P.2, and B.1.1.28 lineages separately. This exploratory analysis revealed similar evolutionary rates within each lineage, but greater root-to-tip distances for P.1 compared to B.1.1.28 (fig. S4), suggesting that the emergence of P.1 was preceded by a period of faster molecular evolution. The B.1.1.7 lineage also exhibits such evolution (13), which was hypothesised to have occurred in a chronically infected or immunocompromised patient (28, 29).

In order to date the emergence of P.1 while accounting for a faster evolutionary rate along its ancestral branch, we used a local molecular clock model (30) with a flexible non-parametric demographic tree prior (31). Using this approach, we estimate the date of the common ancestor of the P.1 lineage to around 6 Nov 2020 (95% Bayesian credible interval, BCI, 9 Oct to 30 Nov 2020). This is approximately one month prior to the resurgence in SARS-CoV-2 confirmed cases in Manaus (Fig. 1A, Fig.2). The P.1 sequences formed a single well-supported group (posterior probability=1.00) that clustered most closely with B.1.1.28 sequences from Manaus (“AM” in Fig. 2), suggesting P.1 emerged there. Further, the local clock model statistically confirms a higher evolutionary rate for the branch immediately ancestral to lineage P.1 compared to that for lineage B.1.1.28 as a whole (Bayes factor, BF=5.25).

Fig. 2.

Visualization of the time-calibrated maximum clade credibility tree reconstruction for B.1.1.28, P.1 and P.2 lineages (n=974) in Brazil.

Terminal branches and tips of Amazonas state are coloured in brown and those from other locations are coloured in green. Nodes with posterior probabilities of <0.5 have been collapsed into polytomies and their range of divergence dates are illustrated as shaded expanses.

Our data suggests multiple introductions of the P.1 lineage from Amazonas to Brazil’s south-eastern states (Fig. 2). We also detected 7 small well-supported clusters of P.2 sequences from Amazonas (2–7, posterior probability=1.00). Virus exchange between Amazonas state and the urban metropolises in southeast Brazil largely follow patterns in national air travel mobility (Fig. 1D, fig. S12).

Infection with P.1 and sample viral loads

We analysed all SARS-CoV-2 RT-qPCR positive results from a laboratory providing testing in Manaus since May 2020 (Fig. 1A, Data S2) with the aim of exploring trends in sample RT-qPCR cycle threshold (Ct) values, which are inversely related to sample virus loads and transmissibility (32). By focusing on data from a single laboratory, we reduce instrument and process variation that can affect Ct measurements.

We analysed a set of RT-qPCR positive cases for which virus genome sequencing and lineage classification had been undertaken (n = 147). Using a logistic function (Fig. 3A) we find that the fraction of samples classified as P.1 increased from 0% to 87% in around 7 weeks (table S1), quantifying the trend shown in Fig 1C. We found a small but statistically significant association between P.1 infection and lower Ct values, for both the E gene (lognormal regression, p = 0.029, n = 128 samples, 61 of which were P.1) and N gene (p = 0.01, n = 129, 65 of which were P.1), with Ct values lowered by 1.43 (0.17–2.60 95% CI) and 1.91 (0.49–3.23) cycles in the P.1 lineage on average, respectively (Fig. 3B).

Fig. 3.

Temporal variation in the proportion of sequenced genomes belonging to P.1, and trends in RT-qPCR Ct values for COVID-19 infections in Manaus.

(A) Logistic function fitting to the proportion of genomes in sequenced infections that have been classified as P.1 (black circles, size indicating number of infections sequenced), divided up into time-periods where the predicted proportion of infections that are due to P.1 is <1/3 (light brown), between 1/3 and 2/3 (green) and greater than 2/3 (grey). For the model fit, darker ribbon represents the 50% credible interval, and lighter ribbon represents the 95% credible interval. For the data points, grey thick line is the 50% exact Binomial confidence interval and the thinner line is the 95% exact Binomial confidence interval. (B) Ct values for genes E and N in a sample of symptomatic cases presenting for testing at a healthcare facility in Manaus, stratified according to the period defined in (A) in which the oropharyngeal and nasal swab collections occurred in. (C) Ct values for genes E and N in a subsample of 184 infections included in (B) that had their genomes sequenced (dataset A).

Using a larger sample of 942 Ct values (including an additional 795 samples for which no lineage information was available) we investigated Ct values across three time periods characterised by increasing P.1 relative abundance. Average Ct-values for both the E and N genes decline through time, as both case numbers and the fraction of P.1 infections increased (Fig. 3C; E gene p = 0.022 and p<0.0001 for comparison of time periods 2 and 3 to period 1; N gene p = 0.025 and p<0.001, respectively).

However, population-level Ct distributions are sensitive to changes in the average time since infection when samples are taken, such that median Ct values can decrease during epidemic growth periods and increase during epidemic decline (33). In an attempt to account for this effect, we assessed the association between P.1 infection and Ct levels whilst controlling for the delay between symptom onset and sample collection. Statistical significance was lost for both data sets (E gene p = 0.15, n = 42, 22 of which were P.1; N gene p = 0.12, n = 42, 22 of which were P.1). Due to this confounding factor we conclude that we cannot yet determine if P.1 infection is associated with increased viral loads (34) or a longer duration of infection (35).

Mathematical modelling of lineage P.1 epidemiological characteristics

We next explored epidemiological scenarios that might explain the recent resurgence of transmission in Manaus (36). To do this, we extend a semi-mechanistic Bayesian model of SARS-CoV-2 transmissibility and mortality (41–43) to include two categories of virus (“P.1” and “non-P.1”) and to allow characteristics such as infection severity, transmissibility and propensity for re-infection to vary between the categories. It also integrates information on the timing of P.1 emergence in Manaus, using our molecular clock results (Fig. 2). The model explicitly incorporates waning of immune protection following infection, parameterized using dynamics from a longitudinal cohort study (37), to explore the competing hypothesis that waning of prior immunity might explain the observed resurgence (36). We use the model to evaluate the statistical support that P.1 possesses altered epidemiological characteristics compared to local non-P.1 lineages. The model is fitted to both COVID-19 mortality data (with a correction for systematic reporting delays (38, 39)) and to the estimated increase through time in the proportion of infections due to P.1 derived from genomic data (Table S1). We assume within-category immunity wanes over time (50% wane within a year, though sensitivity analyses are presented in Table S4) and that cross-immunity (the degree to which previous infection with a virus belonging to one category protects against subsequent infection with the other) is symmetric between categories (see details in Supplementary Materials).

Our results suggest the epidemiological characteristics of P.1 are different to those of previously circulating local SARS-CoV-2 lineages, but also highlight substantial uncertainty in the extent and nature of this difference. Plausible values of transmissibility and cross immunity exist in a limited area but are correlated (Fig. 4A, with the extent of immune evasion defined as 1 minus the inferred cross-immunity). This is expected, because in the model a higher degree of cross-immunity means that greater transmissibility of P.1 is required to generate a second epidemic. Within this plausible region of parameter space, P.1 can be between 1.4–2.2 (50% BCI, with a 96% posterior probability of being >1) times more transmissible than local non-P1 lineages and can evade 25–61% (50% BCI, with a 95% posterior probability of being able to evade at least 12%) of protective immunity elicited by previous infection with non-P.1 lineages (Fig. 4A). The joint-posterior distribution is inconsistent with a combination of high increased transmissibility and low cross-immunity (Fig. 4A). Moreover, our results further show that natural immunity waning alone is unlikely to explain the observed dynamics in Manaus, with support for P.1 possessing altered epidemiological characteristics robust to a range of values assumed for the date of the lineage’s emergence and the rate of natural immunity waning (Tables S3 and S4). We caution that these results are not generalisable to other settings; more detailed and direct data are needed to identify the exact degree and nature of the changes to P.1’s epidemiological characteristics compared to previously circulating lineages.

Fig. 4.

Estimates of P.1’s epidemiological characteristics inferred from a multicategory Bayesian transmission model fitted to data from Manaus, Brazil

(A) Joint posterior distribution of the cross-immunity and transmissibility increase inferred through fitting the model to mortality and genomic data. Grey contours refer to posterior density intervals ranging from the 95% and 50% isoclines. Marginal posterior distributions for each parameter shown along each axis. (B) As for (A) but showing the joint-posterior distribution of cross-immunity and the inferred relative risk of morality in the period following P.1’s emergence compared to the period prior. (C) Daily incidence of COVID-19 mortality. Points show severe acute respiratory mortality records from the SIVEP-Gripe database (71), brown and green ribbons show model fit for COVID-19 mortality incidence, disaggregated by mortality attributable to non-P.1 lineages (brown) and the P.1 lineage (green). (D) Estimate of the proportion of P.1 infections through time in Manaus. Black data points with error bars are the empirical proportion observed in genomically sequenced cases (see Fig. 3A) and green ribbons (dark = 50% BCI, light = 95% BCI) the model fit to the data. (E) Estimated cumulative infection incidence for the P.1 and non-P.1 categories. Black data points with error bars are reversion-corrected estimates of seroprevalence from blood donors in Manaus (2), coloured ribbons are the model predictions of cumulative infection incidence for non-P.1 lineages (brown) and P.1 lineages (green). These points are shown for reference only and were not used to fit the model. (F) Bayesian posterior estimates of trends in reproduction number R for the P.1 and non-P.1 categories.

We estimate that infections are 1.1–1.8 (50% BCI, 81% posterior probability of being >1) times more likely to result in mortality in the period following P.1’s emergence, compared to before, although posterior estimates of this relative risk are also correlated with inferred cross-immunity (Fig. 4B). More broadly, the recent epidemic in Manaus has strained the city’s healthcare system leading to inadequate access to medical care (40). We therefore cannot determine whether the estimated increase in relative mortality risk is due to P.1 infection, stresses on the Manaus healthcare system, or both. Detailed clinical investigations of P.1 infections are needed. The model fits well observed time series data from Manaus on COVID-19 mortality (Fig. 4C) and the relative frequency of P.1 infections (Fig. 4D) and also captures previously-estimated trends in cumulative seropositivity in the city (Fig. 4E). We estimate the reproduction number (Rt) on 07 Feb 2021 to be 0.2 (50% BCI: 0.1–0.3) for non-P.1 and 0.5 (50% BCI: 0.4–0.6) for P.1 (Fig. 4F).

Characterisation and adaptation of a constellation of spike protein mutations

Lineage P.1 contains 10 new amino acid mutations in the virus spike protein (L18F, T20N, P26S, D138Y, R190S, K417T, E484K, N501Y, H655Y, T1027I) compared its immediate ancestor (B.1.1.28). In addition to the abovementioned estimated increase in the rate of molecular evolution during the emergence of P.1, we find, using molecular selection analyses (41), evidence that 9 of these 10 mutations are under diversifying positive selection (Fig. S11).

Three key mutations present in P.1, N501Y, K417T and E484K, are located in the spike protein RBD. The former two interact with human angiotensin-converting enzyme 2 (hACE2) (11), whilst E484K is located in a loop region outside the direct hACE2 interface (fig. S11). Notably, the same three residues are mutated with the B.1.351 variant of concern, and N501Y is also present in the B.1.1.7 lineage. The independent emergence of the same constellation of mutations in geographically-distinct lineages suggests a process of convergent molecular adaptation. Similar to what was observed for SARS-CoV-1 (42–44), mutations in the RBD may increase affinity of the virus for host ACE2 and consequently impact host cell entry and virus transmission. Recent molecular analysis of B.1.351 (45) suggests that the three P.1 RBD mutations may similarly enhance hACE2 engagement, providing a plausible hypothesis for an increase in transmissibility of the P.1 lineage. Moreover, E484K has been associated with reduced antibody neutralisation (6, 46–48) and as RBD-presented epitopes account for ~90% of the neutralising activity of sera from individuals previously infected with SARS-CoV-2 (49), tighter binding of P.1 viruses to hACE2 may further reduce the effectiveness of neutralizing antibodies that are competing with hACE2 to bind the RBD. However, it remains difficult to estimate the contributions of the P.1 RBD mutations to transmissibility and neutralisation that may have led to P.1’s emergence at the population level.

Conclusion

Genomic surveillance and early data sharing by teams worldwide led to the rapid detection and characterisation of P.1 (23), yet such surveillance is still limited in many settings. Existing genomic surveillance is currently inadequate to determine the true international extent of P.1, and this paucity limits the detection of similar variants of concern globally. Studies to evaluate real-world vaccine efficacy in response to P.1 are urgently needed, although we note that neutralisation titres represent only one component of the elicited response to vaccines, and that minimal reduction of neutralisation titres relative to earlier circulating strains is not uncommon (45). Until an equitable allocation and access to effective vaccines is available to all, non-pharmaceutical interventions should continue to play an important role in reducing the emergence of new variants.