Comparative Analysis of Emerging B.1.1.7+E484K SARS-CoV-2 Isolates.

We report the genome of a B.1.1.7+E484K severe acute respiratory syndrome coronavirus 2 from Southeastern Pennsylvania and compare it with all high-coverage B.1.1.7+E484K genomes (n = 235) available. Analyses showed the existence of at least 4 distinct clades of this variant circulating in the United States and the possibility of at least 59 independent acquisitions of the E484K mutation.

During the past 6 months of the pandemic, several variants of concern (VOC), each represented by a constellation of specific mutations thought to enhance viral fitness, have emerged in viral lineages from the United Kingdom (20I/501Y.V1; B.1.1.7), South Africa (20H/501Y.V2; B.1.351), and Brazil (20J/501Y.V3; P.1). These lineages were concerning due to likely increased transmission rates [1-6]. Two of these lineages, B.1.351 and P.1, were of specific concern because they have the mutation E484K, which has been shown to enhance escape from neutralizing antibody inhibition in vitro [7] and may be associated with reduced efficacy of the vaccine [8-11]. In general, viruses from the B.1.1.7 lineage do not have this mutation. However, in February 2021 Public Health England (PHE) published a concerning report of 11 B.1.1.7 genomes that had acquired the E484K spike mutation [12].

Here we report a B.1.1.7 isolate with the E484K spike mutation isolated in Southeastern Pennsylvania (PA). Our laboratory at the Children’s Hospital of Philadelphia performed sequencing on randomly selected isolates collected since January 2021. Figure 1A shows the diversity of 114 randomly sequenced genomes. Lineages B.1.1.7, B.1.429 (California), B.1.526 (New York), and R.1 (international lineage with the E484K mutation) accounted for 69% of the sequenced genomes in March. There was a massive increase in lineage B.1.1.7 from 2% (1/47) in February to 42% in March (15/36). Interestingly, 1 B.1.1.7 isolate from this surveillance, collected on 3/24/2021 from a 52-year-old male, carried the E484K spike mutation that is present in the South African and Brazilian lineages.

Figure 1.

Diversity of SARS-CoV-2 in Philadelphia and global diversity of sequenced B.1.1.7+E484K genomes. A, Stacked area plot showing the diversity of random genomes sequenced by our laboratory at Children’s Hospital of Philadelphia from January, February, and March 2021. Ten lineages that were represented by only 1 genome (B.1.1, B.1.1.106, B.1.1.129, B.1.1.197, B.1.1.281, B.1.1.296, B.1.119, B.1.234, B.1.350, B.1.409) were excluded from the plot. One B.1.526.1 isolate was counted with the parent B.1.526 for easier visualization. B, Stacked area plot showing the number of GISAID genomes (n = 250) that are B.1.1.7 (20I/501Y.V1) and have the E484K spike mutation over time in the United States and globally. C, Diversity of 236 isolates according to GNUVID. Stacked area plot showing relative abundance of circulating CCs for the 236 B.1.1.7+E484K isolates (typed by GNUVID). The bar plot shows that the isolates belong to 7 different CCs. Isolate EPI_ISL_1385215 was not assigned to any of the 7 CCs (CC255). Fourteen isolates were excluded from the plot as they had >5% nucleotides designated “N” in the sequence. Abbreviations: CCs, clonal complexes; GISAID, Global Initiative on Sharing All Influenza Data; SARS-CoV-2, severe acute respiratory syndrome coronavirus 2.

To better understand the relationship between this isolate and publicly available severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) genomes, we compared it with all available B.1.1.7+E484K high-coverage genomes available on the Global Initiative on Sharing All Influenza Data (GISAID; n = 235) [13]. Since the first report by PHE in February, a total of 253 B.1.1.7+E484K genomes have been uploaded to GISAID from England and 14 other countries (Germany, France, Italy, Poland, Sweden, Ireland, Netherlands, Portugal, Wales, Turkey, Slovakia, Austria, Czech Republic, and the United States) [13].

A temporal plot of the number of global B.1.1.7+E484K isolates collected between December 2020 and March 2021 (2-week window) is shown in Figure 1B (as of 04/17/2021). The first isolate of the 60 US isolates available on GISAID [13] was collected on 02/06/2021 from Oregon. Isolates were also reported from 15 other states (New York, North Carolina, Connecticut, Georgia, New Jersey, Maryland, Florida, West Virginia, California, Pennsylvania, Michigan, Texas, Massachusetts, Washington, and Colorado). Of these isolates, 48% were from Florida (n = 17) and New York (n = 12) and 28% were from New Jersey (n = 7), California (n = 4), and Pennsylvania (n = 6). Two isolates were from Oregon (OR), Connecticut (CT), and Maryland (MD), and single isolates are recorded from Georgia (GA), Texas (TX), Massachusetts (MA), Washington (WA), Colorado (CO), West Virginia (WV), Michigan (MI), and North Carolina (NC). The number of US isolates in March (n = 47 including the PA isolates) was nearly 6 times the number of the isolates reported in February. Since this analysis was completed and at the time of this submission (05/26/2021), there are 1400 B.1.1.7+E484K genomes on GISAID [13], which raises the concern that more B.1.1.7+E484K sequences may be emerging even as herd immunity increases by natural immunity and vaccines.

Although all 236 B.1.1.7+E484K genomes were typed as B.1.1.7 using Pangolin [14], a more granular view using our typing tool “GNUVID” [15] shows that they belong to 7 different clonal complexes (CCs; 45062, 46649, 49676, 57630, 58534, 62415, and 67441) (Figure 1C; Supplementary Table 1). In the GNUVID typing system, these correspond to 7 of 10 CCs in the B.1.1.7 lineage. For each of these CCs, representative sequences without the E484K mutation have been circulating since at least November 2020, predating the first E484K in each CC. This raised the possibility that the E484K mutation was acquired independently in each of these CCs in independent events.

To test the hypothesis of multiple acquisitions in a rigorous phylogenetic framework, we queried the GISAID database for closely related genomes for each of the 236 high-quality B.1.1.7+E484K genomes and retrieved 354 closely related B.1.1.7 genomes. Using the phylogeny reconstructed from the combined data set, we performed an ancestral reconstruction using Fitch optimization [16], which revealed 59 de novo acquisitions of the E484K mutation in the B.1.1.7 lineage (Figure 2). In the phylogeny, 6 E484K clades were found in >1 country and 4 E484K clades in >1 US state, suggesting widespread dissemination of these viruses. The genome from PA presented here falls in a well-supported clade of 77 B.1.1.7 isolates, 48 of which have the E484K mutation. Seven of the 48 isolates were from the United States (CT, FL, OR, PA, and NY), 33 from Sweden, 4 from Denmark, 2 from Poland, 1 from the Netherlands, and 1 from Germany (Figure 2). The 9 other B.1.1.7+E484K isolates reported from PA were in a large clade containing US genomes (from CA, CO, FL, MA, MD, MI, NC, NE, NJ, NY, OR, TX, WA, WV). This large well-supported clade also contained isolates from England. Another large clade of 75 B.1.1.7 isolates had 34 isolates that carry the E484K mutation; 30 of them are from the United States (CA, CT, FL, MA, NC, NH, NJ, NY, OR). This clade also contained isolates from Germany, Wales, and Ireland. Together these analyses support sustained transmission of these 3 clades within the United States and also globally.

Figure 2.

SNP-based phylogeny showing the independent acquisitions of E484K in the B.1.1.7 lineage. Maximum likelihood tree of the B.1.1.7+E484K isolates. US isolates are in red. For the CHOP_204 isolate, the alternative allele was called a consensus if its frequency was at least 0.75. The tree was rooted with MN908947.3. The countries of the isolates are shown as a ring. The PA B.1.1.7+E484K isolates are represented with red stars. The red branches represent ancestral reconstruction of the E484K mutation in the B.1.1.7 lineage. Bootstrap values >70 are shown on the branches. Abbreviations: SNP, single nucleotide polymorphism; PA, Pennsylvania.

Phylogenetic and single nucleotide polymorphism (SNP) analysis in the 236 isolates compared with the reference MN908947.3 (Supplementary Figures 1 and 2) [17] showed that the isolate presented here had 12/17 of the B.1.1.7-defining SNPs (Supplementary Table 2), while the other Pennsylvanian isolate in the same clade had 17/17 of the SNPs. It also shared with 9 other US isolates a stop mutation (A28095T) in ORF8 (Supplementary Figure 2).

Here we present a comparative analysis of the first SARS-CoV-2 B.1.1.7 isolates detected in PA that carry the E484K spike mutation, a mutation that could be associated with reduced efficacy of both vaccine-induced and natural immunity. Our analysis shows that multiple lineages of B.1.1.7+E484K are circulating in the United States and globally and that these lineages have acquired the E484K mutation independently, which argues for strong selective pressure for this mutation.

METHODS

A nasopharyngeal swab sample that had residual volume after initial laboratory processing, positive PCR testing for SARS-CoV-2, was obtained for this study. RNA was extracted from nasopharyngeal swab samples using the QIAamp Viral RNA Mini (Qiagen). Whole-genome sequencing was done by the Genomics Core Facility at Drexel University. Briefly, whole-genome sequencing of extracted viral RNA was performed as previously described using the Paragon Genomics CleanPlex SARS-CoV-2 Research and Surveillance NGS Panel [18, 19]. Libraries were quantified using the Qubit dsDNA HS (High Sensitivity) Assay Kit (Invitrogen) with the Qubit Fluorometer (Invitrogen). Library quality was assessed using the Agilent High Sensitivity DNA Kit and the 2100 Bioanalyzer instrument (Agilent). Libraries were then normalized to 5 nM and pooled in equimolar concentrations. The resulting pool was quantified again using the Qubit dsDNA High-Sensitivity (HS) Assay Kit (Invitrogen) and diluted to a final concentration of 4 nM; libraries were denatured and diluted according to Illumina protocols and loaded on the MiSeq at 10 pM. Paired-end and dual-indexed 2×150-bp sequencing was done using MiSeq Reagent Kits, version 3 (300 cycles). Sequences were demultiplexed and basecalls were converted to FASTQ using bcl2fastq2, version 2.20. The FASTQ reads were then processed to consensus sequence, and variants were identified using the ncov2019-artic-nf pipeline (https://github.com/connor-lab/ncov2019-artic-nf). Briefly, the pipeline uses iVar [20] for primer trimming and consensus sequence-making (options: --ivarFreqThreshold 0.75). A bed file for the Paragon kit primers was used in the pipeline.

All 253 SARS-CoV-2 genomes that were assigned to Pango lineage [14] B.1.1.7 and possessing the E484K spike mutation (including the study isolate CHOP_204) were downloaded from GISAID [13] on 04/17/2021. An acknowledgments table of the submitting laboratories providing the SARS-CoV-2 genomes used in this study is in Supplemental Table 3. Seventeen sequences were excluded for lower coverage (>5% Ns; n = 14) and missing collection date (n = 3). All the high-coverage SARS-CoV-2 genomes (n = 236) were assigned a clonal complex using the GNUVID, version 2.2, database [15]. Temporal plots were plotted using a custom script.

To show the relationship among the genomes of the 236 isolates, a maximum likelihood tree was constructed. Briefly, consensus SARS-CoV-2 sequences for the 236 isolates were aligned to MN908947.3 [17] using MAFFT’s FFT-NS-2 algorithm [21] (options: --add --keeplength)). The 5’ and 3’ untranslated regions were masked in the alignment file using a custom script. A maximum likelihood tree using IQ-TREE 2 [22] was then estimated using the GTR+F+I model of nucleotide substitution [23], default heuristic search options, and ultrafast bootstrapping with 1000 replicates [24]. The tree was rooted to MN908947.3. The snipit tool was then used to summarize the SNPs in the 236 isolates relative to MN908947.3 (https://github.com/aineniamh/snipit). To investigate the number of independent acquisitions of the E484K mutation, a maximum likelihood tree was constructed that has both the B.1.1.7 and B.1.1.7+E484K genomes. Briefly, all available genomes were downloaded from GISAID. A Mash [25] database of all available GISAID SARS-CoV-2 genomes was sketched (options: -i -p 32 -k 32 -s 10000). As there are hundreds of thousands of B.1.1.7 genomes available, we identified the 3 closest genomes to each of the 236 B.1.1.7+E484K genomes by sorting the GISAID genomes using the Mash distance (options: -i -p 32 -d 0.00055). This process identified 354 nonredundant B.1.1.7 and B.1.1.7+E484K genomes close to the 236 study genomes. A maximum likelihood tree of the 590 (354 + 236) genomes was then produced as discussed. The trees were visualized in iTOL [26].

Supplementary Data

Supplementary materials are available at Open Forum Infectious Diseases online. Consisting of data provided by the authors to benefit the reader, the posted materials are not copyedited and are the sole responsibility of the authors, so questions or comments should be addressed to the corresponding author.

Click here for additional data file.