Regional Replacement of SARS-CoV-2 Variant Omicron BA.1 with BA.2 as Observed through Wastewater Surveillance.

Greater knowledge of circulating SARS-CoV-2 variants can inform pandemic response, vaccine development, disease epidemiology, and use of monoclonal antibody treatments. We developed custom assays targeting characteristic mutations in SARS-CoV-2 variants Omicron BA.1 and BA.2 and confirmed their sensitivity and specificity in silico and in vitro. We then applied these assays to daily wastewater solid samples from eight publicly owned treatment works in the greater Bay Area of California, United States, over four months to obtain a spatially and temporally intensive data set. We documented regional replacement of BA.1 with BA.2 in agreement with, and ahead of, clinical sequencing data. This study highlights the utility of wastewater surveillance for real-time tracking of SARS-CoV-2 sublineage circulation. The results suggest that concerted efforts to design RT-PCR assays that target variant and variant sublineage characteristic mutations for wide-scale wastewater monitoring implementation will be informative for pandemic response.

Introduction

Knowing when a new SARS-CoV-2 variant arrives locally is of utmost importance to pandemic response to ensure a rapid response to new features of the virus. In late 2021, the Omicron variant emerged in South Africa and quickly spread to Europe and then North America during a period of time during which Delta was the dominant circulating variant. The Omicron variant resulted in a virus more transmissible and less susceptible to therapeutics and vaccines than previous lineages.[1] Given potential implications for public health, timely identification of its introduction into a new community was critical. Variant identification is usually accomplished through sequencing of clinical specimens, but the time between specimen collection and the return of sequencing information to public health officials is approximately 2 weeks. For example, outbreak.info,[2] a website that provides analysis and visualization of SARS-CoV-2 clinical sequences deposited in GISAID,[3,4] stresses data are unreliable and incomplete during the 1–2 week window preceding the query date.

Wastewater represents a biological composite sample from the contributing community. It contains feces, urine, sputum, mucus, blood, vomitus, and any other excretion that goes down a drain into the sewer network. Global investigations have illustrated that SARS-CoV-2 RNA concentrations, whether measured by RT-QPCR or digital RT-PCR or in the liquid or solid phase of wastewater, correlate well with laboratory-confirmed incident COVID-19 cases in the associated sewersheds.[5−8] A number of studies have illustrated that circulating or novel variants can be identified in wastewater using RT-PCR assays that target a characteristic variant mutation[9−11] and by sequencing SARS-CoV-2.[12,13] We previously developed targeted RT-PCR assays for characteristic mutations in Alpha, Mu, Lambda, Delta, and Omicron BA.1 and showed that the concentrations of these mutations, normalized by a pan-SARS-CoV-2 nucleocapsid (N) gene target, are strongly correlated to the fraction of clinical specimens classified as the associated variant[10,14] at two publicly owned treatment works (POTWs). Here, we extend that previous work by developing and testing two new assays, one that targets characteristic mutations in the BA.2 sublineage of Omicron and another that targets a set of mutations present in both BA.1 and BA.2. We then applied the assays to wastewater samples collected daily at eight POTW in the greater Bay Area, California, United States. Results presented here were available publicly within 24 h of sample collection via our website (wbe.stanford.edu), and these data allowed us to observe near real-time regional replacement of BA.1 with BA.2.

Materials and Methods

Characteristic Spike Gene Mutation Assay Development

Digital RT-PCR assays were designed to target the BA.2 characteristic mutation LPPA24S (a 9 bp deletion) and a region of five adjacent SNPs common to BA.1 and BA.2 (S477N, T478K, Q493R, Q498R, Y505H, hereafter S:477–505). The del143–145 assay targeting a mutation present in BA.1 was developed in a similar manner described elsewhere.[10] According to sequences deposited in GISAID and accessed via outbreak.info (accessed 16 April 2022), the LPPA24S mutation is present in 92% of BA.2 genomes, including sublineages (hereafter, BA.2*), and 0.003% of BA.1 genomes, including sublineages (hereafter, BA.1*) globally, whereas the del143–145 mutation is present in 93% BA.1* and 0% BA.2* genomes, globally (Table S1). Assays were developed in silico using Primer3Plus (https://primer3plus.com/) (Table S2).

Primers and probe sequences (Table ) were screened for specificity in silico using NCBI BLAST and then tested in vitro against various respiratory viruses and gRNA from wild-type (WT) SARS-CoV-2 and various SARS-CoV-2 variants (see SI). Viral RNA was used undiluted as a template in digital droplet RT-PCR with mutation primer and probes (see below). The concentration of gRNA used for in vitro specificity testing was approximately 275 copies per well.

Table 1

Primers and Probes for New Mutation Assaysa

Mutations		Sequence
LPPA24S	Forward	GCCACTAGTCTCTAGTCAGTGTG
	Reverse	TGTCAGGGTAATAAACACCACGT
	Probe	CAGAACTCAATCATACACTAATTCTTTCAC (5′ FAM or HEX/ZEN/3′ IBFQ)
S:477–505	Forward	CTATCAGGCCGGTAACAAAC
	Reverse	ACTACTCTGTATGGTTGGTGAC
	Probe	CCTTTACGATCATATRGTTTCCGACCC (5′ FAM or HEX/ZEN/3′ IBFQ

Primers and probes for other assays that are previously published are in Table S5. Information on the fluorescent molecule and quenchers used for the probes are provided in parentheses after their sequence. FAM, 6-fluorescein amidite; HEX, hexachloro-fluorescein; ZEN, a proprietary internal quencher from IDT; IBFQ, Iowa Black FQ.

The sensitivity and specificity of the mutation assays were further tested by diluting target variant gRNA (Table S1) in no (0 copies), low (100 copies), and high (10,000 copies) backgrounds of WT SARS-CoV-2 gRNA.

Wastewater Samples

Eight POTWs serving populations of the Greater Bay Area and Sacramento, California, United States, were included in the study. They serve between 66,000 and 1.5 million people and offer good spatial coverage of the region; further descriptions are provided elsewhere.[15]

Here, 50 mL of settled solids was collected each day from each POTW between 1 January 2022 and 12 or 13 April 2022, during the period of time when Omicron BA.2 emerged in the region. Detailed sample collection methods are provided elsewhere[15] and in the SI. Samples were immediately stored at 4 °C and transported to the laboratory by a courier service where processing began within 6 h of collection. Between 94 and 103 samples from each POTW were included in the analysis (total of 810).

RNA was extracted from 10 replicate aliquots of each dewatered settled solid sample, and PCR inhibitors were removed[15] (see SI). RNA was processed immediately to measure concentrations of the N gene of SARS-CoV-2,[16] the BA.1* characteristic mutation (del143–145),[10] the BA.2* characteristic mutation (LPPA24S), pepper mild mottle virus (PMMoV), and bovine coronavirus (BCoV) recoveries using digital droplet RT-PCR methods.[15,17] The N gene target is present in all variants and is a pan-SARS-CoV-2 gene target. PMMoV is highly abundant in human stool and wastewater globally[18,19] and is used as an internal extraction/process control for the wastewater samples.[20] BCoV was spiked into the samples and used as an additional recovery control; all samples were required to have greater than 10% BCoV recovery. RNA extraction and PCR negative and positive controls were included to ensure no contamination. Extracted RNA samples from one POTW (SJ) were then stored at −80 °C between 1 and 7 days before they were analyzed for S:477–505 and the N gene. Each of the 10 replicate RNA extracts was run in its own well, and the 10 wells were merged for analysis. Wastewater data are available publicly at the Stanford Digital Repository (https://doi.org/10.25740/cf848zx9249); results below are reported as suggested in the EMMI guidelines.[21]

ddRT-PCR

The ddRT-PCR methods applied to wastewater solids to measure PMMoV and BCoV are provided in detail elsewhere.[15] The N gene and del143–145 mutation were multiplexed along with a pan-SARS-CoV-2 S gene assay (not reported herein). The LPPA24S mutation was multiplexed in an assay along with influenza A and RSV gene targets (only LPPA24S reported herein). The S:477–505 mutation was multiplexed with the N gene and a characteristic mutation in Delta (not reported). Each 96-well PCR plates of wastewater samples included PCR and extraction positive and negative controls. Additional details of the PCR methods are provided in the SI.

QuantaSoft Analysis Pro Software (Bio-Rad, version 1.0.596) was used for thresholding. For a sample to be recorded as positive, it had to have at least three positive droplets. For the wastewater solid samples, three positive droplets across 10 merged wells corresponds to a concentration between ∼500–1000 copies per g (cp/g); the range in values is a result of the range in the equivalent mass of dry solids added to the wells.

Concentrations of RNA targets were converted to concentrations per dry weight of solids in units of copies/g dry weight using dimensional analysis. Dry weight of dewatered solids was determined by oven drying.[22] The total error is reported as standard deviations and includes the errors associated with the Poisson distribution and the variability among the 10 replicates.

Variants Present in Regional Clinical Specimens

The 7-day, centered, rolling average proportions of clinical specimens sequenced from California classified as Omicon BA.1* and BA.2* as a function of specimen collection data were acquired from GISAID[3,4] on 21 April 2022.

Statistics

We normalized mutation concentrations by N gene concentrations to represent the fraction of total SARS-CoV-2 RNA (represented by the N gene assay target which is conserved across variants) that comes from the variant, hereafter referred to as relative mutation concentration. We used Kendall’s tau to test for associations between the 5-day trimmed smoothed average relative concentrations and the proportion of clinical specimens assigned to the corresponding variants as the two variables were not normally distributed (Shaprio Wilk test, p < 0.05 for all). The daily measured relative concentration was matched to the 7-day, centered, rolling average fractions of clinical specimens classified as the associated variant obtained from GISAID. All analysis was carried in R Studio version 1.4.1106 using R version 4.0.5.

Results and Discussion

Variant Mutation Assay Specificity

In silico analysis indicated no cross reactivity between LPPA24S (targeting BA.2*) and S:477–505 (targeting BA.1* and BA.2*) assays and sequences deposited in NCBI. When challenged against a respiratory virus panel and various SARS-CoV-2 variant gRNA (including WT, Alpha, Beta, Gamma, Delta, Mu, Lambda, and BA.1; Table S1), no cross reactivity was observed; the assays only amplified gRNA from BA.2 and BA.1 and BA.2, respectively. Positive controls and NTCs run on sample plates were positive and negative, respectively. When target gRNA was diluted into a background of low and high nontarget gRNA from WT SARS-CoV-2, measurements were similar when there was no, low, or high background gRNA (Figure S1). Results indicate the new mutation assays are specific and sensitive. Wolfe et al.[10] showed the del143–145 (targeting BA.1) assay was specific and sensitive.

Variant RNA Concentrations in Wastewater Solids

All positive and negative controls were positive and negative, respectively, indicating assays performed well and without contamination. BCoV recoveries were higher than 10%, and PMMoV concentrations within the expected range for the POTW suggested an efficient and acceptable recovery of RNA during RNA extraction (Figure S2).

Across all eight POTWs, a common pattern in N, del143–145, and LPPA24S concentrations is observed. Between 1 January and late Februrary to mid-March 2022, N and del143–145 concentrations are approximately the same, and both decrease from 106 to 104 cp/g (Figure ). The 2 orders of magnitude decrease is coincident with a decline in laboratory-confirmed incident COVID-19 cases in the sewersheds after the winter Omicron surge in the region (case data not shown). In late February to mid-March until mid-April, both N and LPPA24S concentrations increase from 103 to 105 cp/g, while del143–145 concentrations decrease from 104 to 102 cp/g or nondetect. LPPA24S concentrations are similar to those of N at the POTWs during that same time period.

Figure 1

Concentrations of N, del143–145, LPPA24S, and S:477–505 at eight POTWs in the study (name indicated in corner). Error bars are 68% confidence intervals that include error from variation among 10 replicate wells and Poisson error (“total error” as reported by the instrument). If error bars cannot be seen, they are smaller than the symbol. Asterisks near the x-axis indicate samples that are nondetects for the associated target as indicated by the color.

The ratios del143–145/N and LPPA24S/N describe the concentrations of the mutations in relation to a conserved gene target. As a mutation’s relative concentration approaches 1, it suggests that the variant associated with the mutation may be dominant in the wastewater. The Del143–145 relative concentration was a maximum and close to 1 at most POTWs at the beginning of the time series (1 January 2022) and then fell to nearly 0 by mid-April 2022 (Figure ). At the same time, the LPPA24S relative concentration increased from 0 when it was first measured to near 1 at all POTWs. Depending on POTW, the date on which the LPPA24S relative concentration surpassed that of del143–145 ranged from 26 February to 19 March (Table S3). This date was earliest (26 February) at Ocean, located in San Francisco, and latest at Gil and Sac; the date is positively correlated with the distance from Ocean (r = 0.9, p < 0.05, n = 8), suggesting a spatial pattern in BA.2* emergence. Regional replacement of BA.1* with BA.2* appears to be complete at all POTWs by the end of the times series as the del143–145 relative concentrations are close to 0, and LPPA24S relative concentrations approach 1. At times, relative concentrations may exceed 1; this may be due to variability in the performance of the digital RT-PCR assays. It is important to note that the relative concentrations do not incorporate measurement error of the numerator or denominator.

Figure 2

Relative concentrations of del143–145 and LPPA24S (their concentrations normalized by the N gene) at eight POTWs in the study (name indicated in corner). Markers are raw ratio values. Lines represent 5-day trimmed moving averages.

At SJ, we measured the concentration of an additional set of adjacent mutations present in both BA.1 and BA.2 (S:477–505) (Figure ). S:477–505 concentrations agree well with the N concentrations during a period of time when BA.1* was dominant (prior to mid-March) and when BA.2* was dominant (thereafter).

The del143–145 and LPPA24S relative concentrations at each POTW are positively and significantly correlated with the fraction of clinical cases in California assigned as BA.1* and BA.2*, respectively (Table S4). The tau between relative del143–145 concentrations and the proportion of cases assigned as BA.1* ranged from 0.51 at Gil to 0.81 at SAC (median tau = 0.68, all p < 10–7, n = 94–103). The tau between relative LPPA24S concentrations and the proportion of cases assigned as BA.2* ranged from 0.67 at Gil to 0.90 at SAC (median tau = 0.75, all p < 10–5, n = 58–62). Tau values were similar when calculated using raw relative concentrations (see SI). State-level clinical data suggest that the proportion of cases caused by BA.2* exceeds that of BA.1* after 21 March (Figure ), later than any wastewater suggests this occurred regionally. Ideally, these analyses should be performed using data on clinical variants in the sewersheds of the eight POTWs. However, such data are not readily available, so we had to rely on state-level data.

Figure 3

Percentage of clinical cases in California classified as indicated BA.1* and BA.2* as a function of specimen collection date; data represent 7-day moving average as acquired from GISAID. Data available through 13 April when downloaded from GISAID on 21 April 2022. Gray area represents the time period for which data were incomplete on 14 April 2022, when wastewater measurements for the last date in our data series (12 or 13 April 2022) were available.

State-level clinical sample sequencing data on variant occurrence are typically 1–2 weeks delayed owing to the time required for laboratory processing and reporting. In addition, there can be large uncertainties in the proportion estimates preceding the data for which the most recent data are available. This is illustrated in Figure where we highlight the approximate time period for which data were incomplete on 14 April 2022, when wastewater measurements for the last date in our data series (12 or 13 April 2022) were available publicly on our website (wbe.stanford.edu). During this same time period, wastewater data suggest BA.2* is near 100% and BA.1* near 0%, which a week later (when data were downloaded from GISAID on 21 April 2022) were confirmed by clinical sequencing data. This highlights an important advantage of using wastewater surveillance for variant tracking using targeted assays: wastewater measurements obtained using these assays are available within 24 h of receiving samples at the laboratory providing information ahead of clinical sequencing data.

Real-time estimates of variant abundance with high spatial resolution, such as those from wastewater surveillance, could allow clinicians and regulators to identify the most effective treatment for a particular community. The monoclonal antibody therapy sotrovimab provides a recent and practical example. Sotrovimab was granted emergency use authorization by the US FDA for treatment of infection by variants including BA.1, but not BA.2. As BA.2 took over, its use was discontinued by multistate HHS region blocks, informed by clinical sequencing data. While sotrovimab was allowed for use in California until 20 March 2022, our data suggest BA.2 became the dominant variant in the Bay Area weeks prior but before BA.2 had taken over all of California and HHS region 9. Real-time wastewater variant estimates could have allowed the drug’s use to be targeted spatially and temporally to maximize effectiveness.