Skip to main content
Advertisement
  • Loading metrics

Norovirus GII wastewater monitoring for epidemiological surveillance

Abstract

While the Centers for Disease Control and Prevention coordinates several outbreak and clinical surveillance systems for norovirus, norovirus is strongly under-reported due to individuals not seeking care or not being tested. As a result, norovirus surveillance using case reports and syndromic detection often lags rather than leads outbreaks. Digital epidemiology sources such as search term data may be more immediate, but can be affected by behavior and media patterns. Wastewater monitoring can potentially provide a comprehensive and consistent data stream that can help to triangulate across these different data sets. To assess the timeliness of norovirus wastewater testing compared with syndromic, outbreak and search term trend data for norovirus, we quantified human norovirus GII in composite influent samples from 5 wastewater treatment plants (WWTPs) using reverse transcription-digital droplet PCR and correlated wastewater levels to syndromic, outbreak, and search term trend data. Wastewater human norovirus (HuNoV) GII RNA levels were comparable across all WWTPs after fecal content normalization using Pepper mild mottle virus (PMMoV). HuNoV GII wastewater values typically led syndromic, outbreak, and search term trend data. The best correlations between data sources were observed when the wastewater sewershed population had high overlap with the population included by other monitoring methods. The increased specificity and earlier detection of HuNoV GII using wastewater compared to other data, and the ability to make this data available to healthcare, public health, and the public in a timely manner, suggests that wastewater measurements of HuNoV GII will enhance existing public health surveillance efforts of norovirus.

Introduction

SARS-CoV-2 demonstrated the utility of wastewater monitoring of pathogens, and catalyzed the infrastructure development needed to monitor other pathogens of public health importance. One such pathogen is norovirus, which causes more than 90% of epidemic non-bacterial gastroenteritis outbreaks [1]. Norovirus is estimated to result in 900 deaths, 109,000 hospitalizations, 465,000 emergency department visits, and 2.3 million clinic visits [2], annually in the United States, causing substantial health and economic burdens [3].

Norovirus, a genera in the Caliciviridae family, is a nonenveloped virus with a positive-sense RNA genome. There are ten genogroups of Norovirus, GI–GX, and 48 genotypes [4]. The GII genogroup, specifically variants of the GII.4 genotype, are the most common cause of norovirus disease worldwide [3]. The most common route of transmission is person-to-person, followed by contaminated food or water, and contact with contaminated fomites [5]. An estimated 30% of all norovirus infections and over 40% of GII.4 infections are asymptomatic [6]. Reinfection is common due to relatively short-term immunity and the diverse genogroups [7,8]. Norovirus is a “winter” pathogen, in the northern hemisphere most outbreaks occur between November and April [9].

Human norovirus (HuNoV) surveillance practices vary greatly across the US, and there is no requirement for local, territorial, or state agencies to report individual norovirus cases to the national system. Public health officials often rely on syndromic data, such as school and emergency department gastrointestinal illness reports, as early indicators of norovirus outbreaks. This syndromic data is not specific to norovirus and is often only easily available to public health officials in a timely and local manner. Health departments are encouraged to report all waterborne, foodborne, and enteric disease outbreaks, which would include norovirus outbreaks, to the National Reporting System (NORS) and norovirus outbreaks to Calicinet [10]. The Centers for Disease Control and Prevention (CDC) and health departments from fourteen states, including Michigan, participate in the Norovirus Sentinel Testing and Tracking (NoroSTAT) network. This data is specific to norovirus and is made publicly available; however, the lag between norovirus testing and reporting can be as long as a few weeks–unacceptably long for a highly transmissible virus. Furthermore, data aggregated across participating states has minimal value for informing communities and much of norovirus transmission goes undetected by conventional surveillance systems due to asymptomatic cases and cases that do not require clinical care. Research on alternatives to traditional surveillance includes wastewater monitoring and digital epidemiology [11].

Wastewater monitoring for HuNoV has the potential to provide more local, early-warning information to inform public health decision-making—potentially prior to clinically detected outbreaks. Additionally, wastewater data can easily be made publicly available. Unlike other surveillance data, it is not as biased by human care-seeking behavior and clinical testing. However, rather than replace conventional epidemiological monitoring methods, wastewater data provides an alternative data source to triangulate with existing imperfect clinical data streams.

HuNoV is detectable in wastewater treatment plant (WWTP) influent [1214] and studies have reported positive correlations between wastewater levels and conventional surveillance data, including confirmed HuNoV cases [1517], gastrointestinal (GI) illness cases [18,19], and confirmed hospital cases [19,20]. A study in Japan found stronger correlations between wastewater levels and regional rather than national cases [18]. A recent study in the US showed population weighted average HuNoV GII RNA levels across 145 WWTPs positively correlated to national NoroSTAT clinical data [21]. Wastewater data led conventional surveillance data in some reports and not in others [16,18,20]. These discrepancies could be due to a number of reasons, including that the regions covered by the sewershed did not always correspond to the regions covered by the conventional surveillance data, differences between the study populations, the different types of surveillance data used between the studies (e.g., hospital HuNoV cases, confirmed GI cases, syndromic data, etc.), and the low temporal resolution of wastewater and conventional surveillance data used in some of the studies. This suggests a need to examine corresponding wastewater and epidemiological surveillance data for multiple populations at high resolution and concurrently against multiple sources of conventional norovirus surveillance data.

Another new approach to understanding trends for unreported illnesses is digital epidemiology—the use of search data, social media, mobile phone networks, and other digital data not generated for a public health purpose to understand epidemiological patterns [22]. Digital epidemiology has increasingly been used to understand seasonal and other temporal trends of disease across locations using search term data [23,24], suggesting it could also provide a useful data source for early warning of norovirus increases. Google searches have been correlated with norovirus case data [11]; however, this data can be biased based on terminology, behavior, media coverage and access and often lacks spatial resolution down at the community level. This type of digital epidemiological data has not previously been compared with wastewater data, opening up the potential to build a comprehensive early warning system that bridges across multiple clinical, digital, and wastewater data sources.

In this study, we explore cross-correlations between HuNoV GII wastewater data and syndromic, outbreak, and Google search term data at multiple spatial scales and locations over a one year period, enabling us to explore the potential for triangulation across clinical, outbreak, digital, and wastewater data to understand human norovirus patterns. We measured HuNoV GII levels at high temporal resolution in five WWTPs with different sizes and urban/rural composition, analyzing between two and seven samples per WWTP per week when HuNoV GII levels were elevated. We focused on HuNoV GII testing because 392 of the 440 US confirmed norovirus outbreaks recorded by Calicinet (or 89.1%) in 2021–2022 were HuNoV GII, and HuNoV GII has been shown to be more abundant than HuNoV GI in wastewater in North America [12]. Overall, this study assessed the ability of wastewater monitoring of HuNoV GII to provide added value to public health surveillance in combination with other traditional epidemiological surveillance methods.

Methods

The Michigan SARS-CoV-2 Epidemiology–Wastewater Evaluation and Reporting (SEWER) Network includes local partnerships between wastewater utilities, health departments, tribal communities, universities, and laboratories that are collecting and analyzing wastewater for SARS-CoV-2 across the state of Michigan. All SEWER network laboratories use influent samples and approximately the same RNA purification methods. This infrastructure allowed us to monitor samples for additional RNA viruses.

Sample collection procedure

Samples were provided by five WWTPs in southeast Michigan (Ann Arbor, Flint, Jackson, Tecumseh, and Ypsilanti; Table 1) based on agreements established in June 2021. WWTP personnel collected daily influent samples between 7/18/2021 and 7/14/2022, except Tecumseh where sample collection began 1/12/2022. Samples were collected by 24-hour composite samplers kept at 4°C, and 50 ml aliquots were delivered biweekly by courier on ice to the University of Michigan. Samples were stored at 4°C and processed within 120 hours of collection, with limited degradation of target RNAs expected over this time [2527]. Information about each WWTP was obtained from the state of Michigan SWEEP website: https://www.michigan.gov/coronavirus/stats/wastewater-surveillance/dashboard/sentinel-wastewater-epidemiology-evaluation-project-sweepwebsite.

Sample processing

We concentrated viruses in 40 ml samples 8–40 fold using PEG precipitation as previously described with minor modifications to scale down the starting material [28]. Briefly, 50 ml of influent sample was pasteurized for 30 minutes at 56°C. PEG 8000 (Fisher Scientific #BP2331) and 5 M NaCl solution (Sigma-Aldrich #S6546) were added to 40 ml of pasteurized influent samples to final concentrations of 8% (w/vol) and 0.2 M respectively. Samples were gently mixed and incubated at 4ºC for at least 8 hours. Samples were spun down at 4700 x g for 45 minutes at 4ºC. The supernatant was removed, leaving 1–5 mls of concentrated precipitate, with a bias towards lower volumes for samples with less solids to increase likelihood of detection in samples with low virus levels. RNA extraction was performed with 200 ul of sample concentrate using the QIAmp Viral RNA Mini Kit (Qiagen Sciences, MD) with an elution volume of 80 ul of water. A single extraction was performed for each sample collected. Bovine coronavirus (BCoV) (Bovilis Coronavirus Vaccine, Merck Animal Health, NJ) was added to samples prior to RNA extraction reactions as a recovery control. Negative extraction controls and positive BCoV and HuNoV GII extraction controls were prepared. Non-infectious intact HuNoV GI and GII particles (Cat# NATNOV-6MC, Zeptometrix, Buffalo, NY) were used to confirm the efficacy of norovirus RNA extraction and the specificity of the HuNoV GII primer/probe set in ddPCR reactions. PEG precipitation and RNA extraction using the QIAmp Viral RNA kit have previously been shown to be effective methods for norovirus RNA purification [29]. Pepper mild mottle virus (PMMoV) was quantified as a measure of fecal content [30].

RT-ddPCR

Assay details for RT-ddPCR are provided as per MIQE guidelines [31]. Reverse transcription—digital droplet PCR (RT-ddPCR) analysis of HuNoV GII was performed on new RNA samples using multiplexing (with SARS N1 and N2) or triplexing (with BCoV and PMMoV) from March 25, 2022 to July 14, 2022. RNA from samples July 18, 2021 to March 25, 2022 were freeze-thawed once prior to ddPCR, and had been stored at -80°C for 1–8 months. A total of 643 samples were tested, only 3 samples had no detectable HuNoV GII. Gene copies were quantified through one-step RT-ddPCR (n = 3) using the One-step RT-ddPCR Advanced kit for Probes (catalog #1864021, Bio-Rad, CA) with the method described by Flood et al with some modifications [28]. RT-ddPCR reactions were run at 50C for 60 minutes, 95C for 10 minutes, 40 cycles of 95C for 30 seconds and 56C for 1 minute, 98C for 10 minutes, and held at 4C. HuNoV GII-specific primers and probes (synthesized by Integrated DNA Technologies, IA) that target a 97 bp region of the ORF1-ORF2 region and were previously tested for ddPCR were used [32,33]. When HuNoV GII levels were low, 5 μl undiluted RNA was analyzed by multiplexing with the SARS-CoV-2 nucleocapsid 1 (N1) and nucleocapsid 2 (N2) gene targets designed by the US Centers for Disease Control and Prevention (CDC) [34] using the 5’-FAM tagged HuNoV GII probe at half the concentration of the N1 and N2 probes (125 nM for the HuNoV GII probe, 250 nM for SARS-CoV-2 probes). When HuNoV GII levels were high, samples were diluted 1:100 and 5 μl of the diluted RNA was analyzed by triplexing with BCoV and PMMoV targets using both 5’-FAM and 5’-HEX HuNoV GII probe at half the concentration of the BCoV and PMMoV probes (125 nM for both HuNoV GII probes, 250 nM for BCoV and PMMoV probes). Primer concentrations were all 900 nM. All sample and control reactions were run in triplicate. Droplet analysis was performed on the Bio-Rad QX200 droplet digital PCR systems (Bio-Rad, CA) with explanations of thresholding provided in S1 Text. The results from replicate wells were merged. Controls were run along with all PCR reactions and included non-template controls, extraction controls, and positive controls. The N1/N2 primer/probe stocks and positive PCR controls were provided by MSU Rose lab as part of the SEWER project. PMMoV detection methods had previously been established in the laboratory [35] with a published gene block control [30]. A 390 bp G-block (IDT) DNA, specific to the ORF1-2 region of Genebank sequence MT474038.1, was used as a PCR control for HuNoV GII (Table A in S1 Text). To test both the dynamic range of our target and the presence of inhibition in our samples, a dilution series of HuNoV GII DNA from 2,500 to 25 gene copies was added to both nuclease free water and wastewater sample RNA extract prior to performing ddPCR. The same sample RNA extract was tested without the HuNoV GII DNA spiked in to determine the background HuNoV GII RNA levels in the sample. We saw a clear relationship between the amount of HuNoV GII DNA spiked in the extract and the quantities detected in the blank and sample extracts. This indicated the absence of inhibition in the wastewater RNA extracts. Plates containing negative and no template control wells that had greater than 3 droplets were rerun or re-extracted, consistent with the protocol used by the SEWER network for SARS-CoV-2 ddPCR. Samples with less than 30% of the control reaction BCoV values were flagged for further analysis to determine if RNA degradation had occurred and re-extracted from PEG concentrates when appropriate. All primers, probes, and the HuNoV GII DNA control used are listed in Table A in S1 Text [30,32,34,36,37]. The limit of detection was established by serial dilution of the HuNoV GII DNA control and taking into account the limit of the blank (up to 3 positive droplets). Due to differences in sample concentration volumes after PEG precipitation (from 1–5 mls) the limit of detection varied from 6 x 103 to 2.4 x 104 gc/L.

Syndromic, outbreak, and internet search term data

School-reported GI illnesses and emergency department GI-related visit data were provided by local public health departments and are covered through a data use agreement between the Michigan Department of Health and Human Services (MDHHS) and the Wigginton and Eisenberg labs at the University of Michigan (DUA#: 23-UFA01896 with MDHHS). This project was determined to not be regulated by the UM institutional review board (IRB# HUM00218874). These syndromic data included only reports of symptoms associated with norovirus and no tested or confirmed norovirus cases. The collection and analysis methods used in this study complied with the terms and conditions for the source of the data. Per the CDC National Center for Health Statistics standard, to prevent disclosure of individual information, only data with > 10 cases was included. Weekly total norovirus outbreak values from September 25, 2021 through June 12, 2022 were obtained from the CDC NoroSTAT website (www.cdc.gov/norovirus/reporting/norostat) which compiles and reports outbreak data from 14 states including Michigan, representing about 25% of the US population.

Google Trends records how often a term is searched for in a given region relative to its total search volume with the normalized value in a range from 0 to 100. We obtained time-series data for the search terms “norovirus”, “gastroenteritis”, “stomach flu” for the Metropolitan Detroit region and the state of Michigan from 1/1/2021-11/29/2022. Data were joined by wastewater sample collection date and search term data week. Search term values were rolled forward to apply to all dates in a given week for matching purposes.

Statistical analyses

We conducted statistical analyses with Graphpad Prism 9.4.1 except where otherwise noted. We determined the median and interquartile range of log-transformed weekly average HuNoV GII wastewater concentrations (gc/L) and the log-transformed weekly average HuNoV GII/PMMoV wastewater concentration ratios from each WWTP, from all WWTPs combined, and by season.

To ensure a consistent set of methods and to make full use of the numerical values of the data, we used a Pearson correlation for all analyses, although we also evaluated a subset of variables with both Pearson and Spearman correlations and found similar results (files included in our public repository). We calculated Pearson correlations and cross correlations between wastewater HuNoV GII values and the number of school-reported GI illnesses, emergency department GI visits, and weekly total outbreaks from NoroSTAT. Negative cross correlation lag times indicate clinical data leads wastewater data and positive values indicate wastewater values lead (See S1 Text for details).

Pearson correlations and cross correlations were also computed between HuNoV GII against Detroit and Michigan-wide search term trends using R version 4.0.3 (2020-10-10)—"Bunny-Wunnies Freak Out"; Copyright (C) 2020 The R Foundation for Statistical Computing; Platform: x86_64-w64-mingw32/x64 (64-bit).

Result and discussion

HuNoV GII RNA levels and trends for the five individual WWTPs and all WWTPs combined

Wastewater HuNoV RNA levels were highest in the winter and spring, and lowest in the summer and fall (Fig 1). Including all 5 WWTPs, the daily values ranged from 6.0 x 103 gc/L to 9.7 x 107 gc/L (3.8 to 8.0 log10 gc/L) and the median weekly value was 6.0 log10 gc/L, with an interquartile range (IQR) of 5.5–6.5 log10 gc/L (Fig A in S1 Text). Although Tecumseh (TM) had the highest overall values, TM data was only collected between January-July 2022. If TM was excluded the median value is 6.0 log10 gc/L (IQR = 5.5–6.4). These values are higher than what was seen in a previous meta-analysis at 4.7 log10 gc/L (SD 1) for WWTPs in North America [12], but it is unclear if this is a product of the specific year analyzed, method used for analysis, or location.

thumbnail
Fig 1. HuNoV GII wastewater levels at 5 wastewater treatment plants in 2021–2022.

A. Comparison of PMMoV-normalized HuNoV GII levels from 5 WWTPs in Michigan from 2021–2022. HuNoV was quantified in influent samples using ddPCR at least weekly and PMMoV normalized values were plotted over one year. The only exception is TM, where sample collection began later, in January 2022. Lines connect all values measured for each WWTP. B. Seasonal HuNoV GII values in log10 (median, IQR) for each individual WWTP are presented using box-and-whisker plots for all 4 seasons. Note, no fall values for TM were obtained. The definition of the seasons is meteorological, beginning on the 1st day of the months of the equinoxes or solstices.

https://doi.org/10.1371/journal.pwat.0000198.g001

To account for variations in fecal content, we normalized HuNoV GII values to PMMoV values. The general patterns over the year were similar using either HuNoV gc/L concentrations or PMMoV-normalized concentrations (Fig 1 and Fig A in S1 Text). The PMMoV-normalized average weekly values correlated strongly with the average weekly HuNoV gc/L values for samples from each WWTP, with correlation coefficients between 0.72 (AA WWTP) to 0.83 (TM WWTP) (Fig B in S1 Text). In the case of FL, the median gc/L HuNoV GII value was lower than the other WWTPs, and the IQR did not overlap all of the other WWTPs (Fig A in S1 Text). Normalization of FL HuNoV GII values with PMMoV resulted in FL having an overlapping IQR with the other WWTPs, this is likely because FL has a combined sewershed with a broader range of wastewater strength (Fig A in S1 Text). We focus on PMMoV-normalized HuNoV data for the remainder of the manuscript but include the HuNoV GII values in gc/L in the S1 Text to enable direct comparison with previously published data.

When data from the five WWTPs were combined, the highest HuNoV RNA values were obtained in the spring (average values of 0.027 HuNoV GII/PMMoV or 6.6 log10 gc/L) and the lowest in the fall (0.0012 HuNoV GII/PMMoV or 5.4 log10 gc/L); (Fig 1B, Figs C and D in S1 Text). In general, fall and summer HuNoV RNA levels were lower and the winter and spring HuNoV RNA levels were higher for all five WWTPs (Fig 1B, Figs C and D in S1 Text). The HuNoV RNA concentration variations observed within a plant also exhibited seasonal trends. Specifically, the interquartile ranges (IQRs) observed in winter and spring for a specific plant were larger than the IQRs observed in summer and fall for most WWTPs (Fig 1B). The largest variation in HuNoV GII concentration was detected at our smallest WWTP, which is consistent with previous work reporting wastewater plants with smaller sewershed populations having a higher sensitivity to infection clusters [38].

Overall, WWTPs varied by size of sewershed, population and type of sewer system (combined versus separate) (Table 1), but after normalization using PMMoV as a fecal indicator, levels were comparable across all WWTPs. Seasonal trends were consistent with previous reports and showed strong seasonal differences across two orders of magnitude, perhaps partially explaining previous observations of high variation in influent HuNoV RNA levels across studies and over time [12,14], as previous studies typically used lower temporal resolution and shorter time series.

Comparison of HuNoV GII wastewater levels to school GI illnesses and GI-related emergency room visits

To understand the added value of wastewater data to conventional syndromic data, we compared wastewater data with both school-reported gastrointestinal (GI) illnesses and GI-related emergency department (ED) data used by public health partners. From our public health partners we obtained school-reported GI illness and GI-related ED data for Washtenaw and Lenawee counties which were correlated with HuNoV GII wastewater values from the appropriate WWTPs (AA and YC for Washtenaw and TM for Lenawee). Zip code level school-reported GI illness and GI-related ED data was obtained only for correlations with TM wastewater values. State level GI-related ED data was correlated with wastewater HuNoV GII values for all 5 WWTPs. These analyses provided an opportunity to study both how wastewater correlates and also how it leads or lags syndromic data collected at various spatial resolutions.

Wastewater HuNoV GII RNA concentrations lead reported GI illnesses.

Wastewater values and school GI illnesses followed similar trends for TM (Fig 2A), AA and YC WWTPs (Fig 2B). For TM, the school GI numbers are for the town of Tecumseh, whereas for AA and YC, the GI numbers are for all of Washtenaw County. Of note, the AA WWTP sewershed is entirely included within Washtenaw County, whereas only one-third of the YC WWTP sewershed population is within Washtenaw County. Temporal cross-correlations of weekly total school reported GI illnesses or attendance normalized GI illnesses (details in S1 Text), and weekly average HuNoV GII/PMMoV wastewater values highlight the leading signal of wastewater levels compared to cases. Specifically, the highest correlations were observed at +2 weeks for TM (r = 0.82), +3 weeks for AA (r = 0.67) and +3 weeks for YC (r = 0.74; Table 2). The same general correlation trends were observed for each WWTP when HuNoV GII concentrations (gc/L) were used in the correlation analysis instead of PMMoV-normalized values (Fig E in S1 Text). Additionally, because this study uses retrospective data, this lead trend may even be underestimated from a public data availability perspective, as many state and national monitoring systems may only update data on a weekly (or even less frequent) cadence.

thumbnail
Fig 2. Correlation of TM, AA, and YC HuNoV GII wastewater values with syndromic data.

A. TM WWTP HuNoV GII/PMMoV values along with weekly school reported gastrointestinal illnesses at TM schools. B. AA and YC wastewater levels as HuNoV GII/PMMoV and county-level school reported GI illnesses normalized to enrollment. The GI reports per 1000 people were calculated to normalize for large variations in school attendance. C. Correlation between the weekly total number of zip code level GI-related ED visits and TM PMMoV-normalized HuNoV GII wastewater values. D. Correlations between the weekly total number of county-level GI-related ED visits and AA and YC PMMoV-normalized HuNoV GII wastewater values.

https://doi.org/10.1371/journal.pwat.0000198.g002

thumbnail
Table 2. Cross correlations of HuNoV GII wastewater values with conventional monitoring data.

https://doi.org/10.1371/journal.pwat.0000198.t002

Wastewater HuNoV GII RNA concentrations led GI-related emergency department visits.

The HuNoV GII wastewater levels and GI-related ED visits showed similar trends (Fig 2C and 2D), with values rising in early 2022 and peaking in the spring months. When comparing wastewater concentrations with ED visits in the TM zip code and county, the strongest correlation occurred for the same day (r = 0.68 for zip code; r = 0.75 for county; Table 2). Similar results (±1 week) were obtained when we applied the percent of all ED visits that were GI-related rather than total GI-related ED visit numbers. As we reduce the spatial resolution of the GI-related ED visits from zip code- and county-level data to state-level ED data, the cross-correlations have multiple time frames with high correlations, with the highest correlation at +2 weeks (r = 0.88; Table 2).

For AA and YC, we also observed a significant correlation between wastewater levels and total weekly GI-related ED visits in the respective county (Fig 2D). The highest correlation values were obtained at +2 weeks for both WWTPs (AA R = 0.58; YC R = 0.81) (Table 2). When we reduced the public health data resolution to the state level ED GI-related visits, the timing of the highest correlations varied. For AA, the highest correlation was at +3 weeks (R = 0.67) and for YC the highest resolution was at +1 week (R = 0.85; Table 2). Most WWTPs had extended time periods over which a significant correlation was seen when correlated to the state level data, e.g., YC HuNoV GII wastewater values had significant correlation with state syndromic data from -2 to +3 weeks (R = 0.74 to 0.85; Table 2) and FL wastewater values showed a similar pattern with high correlations from 0 to +3 weeks (R = 0.74 to 0.75; Table 2). Similar cross-correlation results were observed for ED data when wastewater HuNoV GII values in gc/L were used instead of HuNoV GII/PMMoV (Fig F in S1 Text). In some cases, correlations were poorer with gc/L; e.g., AA data resulted in no significant correlations.

Overall, this data of wastewater concentrations versus GI-related ED visits demonstrates that wastewater HuNoV GII signals lead GI-related ED signals. Previously published studies showed wastewater norovirus data leading conventional surveillance data in some reports and not in others [16,18,20]. Likely one of the variables affecting these results are differences in overlap between the sewershed population and the norovirus syndromic or case data population. Data sets that aggregate large populations (e.g., the state ED GI visits) incorporate many localized outbreaks that occur seasonally but not at precisely the same times. Consequently, the data sets exhibit broader, more hill-like peaks. By contrast, data sets representing smaller population sizes may be more likely to discern individual or clusters of outbreaks, which for norovirus can be quite rapid/explosive [39] leading to sharper, more ‘spiky’ increases. Our results suggest that smaller populations and closer overlap between the wastewater and syndromic or case populations resulted in obvious peak correlations with wastewater leading the syndromic and case data.

Comparison of HuNoV GII wastewater levels to publicly available data: NoroSTAT outbreaks

In addition to the school GI illness data and GI-related ED visit data, which has high spatial resolution but is difficult to access as it is often only available to local/state public health, we also compared HuNoV GII data to the publicly accessible NoroSTAT outbreak data. NoroSTAT data has relatively poor spatial resolution compared to school GI illnesses and ED-visits; it encompasses norovirus outbreaks for 14 states covering a population of over 80 million individuals, whereas the school GI illnesses and ED visits were for the zip code and county level. The temporal patterns for the 2021–2022 norovirus outbreaks reported by NoroSTAT and HuNoV GII wastewater values were generally similar (Fig 3A). Wastewater data from the 5 WWTPs correlated well with NoroSTAT data (R = 0.50 to 0.94) for each WWTP from -3 to +2 weeks when comparing the timing of the outbreaks to the HuNoV GII wastewater levels (Fig 3B and Table 2). The highest correlation occurred from -2 to +1 weeks, depending on the WWTP. Similar to the state syndromic data, the NoroSTAT data has low spatial resolution; this likely explains the range of positive and negative lag times with high correlations observed compared to community wastewater data.

thumbnail
Fig 3. Correlation of HuNoV GII wastewater values with NoroSTAT reported outbreaks.

A. Time course from August 2021 to June 2022 showing weekly average HuNoV GII/PMMoV values from 4 WWTPs and weekly norovirus outbreaks reported by NoroSTAT. TM WWTP data is only from January–June 2022. B. Cross correlations between NoroSTAT outbreak numbers and weekly average HuNoV GII/PMMoV values. Graph showing Pearson’s correlation coefficients (r) and probability values (p) determined by comparing weekly average HuNoV GII/PMMoV wastewater values from all 5 WWTPs to the weekly total outbreaks reported by NoroSTAT. Cross correlations were tested for a lead (-4 to -1 weeks), same timing (0), or a lag (1 to 3 weeks) in outbreak reports compared to wastewater values. Note: * = p<0.05, ** = p<0.01, *** = p<0.001.

https://doi.org/10.1371/journal.pwat.0000198.g003

Search term data

Google trends provides public access to search term volume for different regions, and from a digital epidemiology perspective, a useful tool for syndromic surveillance. We compared the number of weekly internet searches for terms related to GI illnesses in Michigan and the Detroit metropolitan (metro) area to HuNoV GII wastewater levels at our five WWTPs. An increase in “norovirus” and “stomach flu” searches in early 2022 appeared to coincide with increasing wastewater HuNoV GII levels (Fig 4A shows YC as an example). Further, “stomach flu” in the Detroit metro area and Michigan had positive associations with HuNoV GII wastewater values at all 5 WWTPs, ranging from 0.53 to 0.71 for the Detroit metro area and 0.66 to 0.82 for the state of Michigan (Fig 4B). Four of the 5 plants showed similarly high correlations for the search term “norovirus”. The smallest WWTP, TM, showed lower values (0.38 for the Detroit metro area, 0.46 for MI), potentially due to the small population and earlier timing of local outbreaks (Fig 4B). The highest correlation values were seen for the largest WWTP, YC, which also has its sewershed population in the Detroit metro area included in the google trend search (Fig G in S1 Text). Cross correlations showed variations in timing for HuNoV GII wastewater concentrations compared to search term trends, with the majority, but not all, showing wastewater values leading search term trend data (Figs H and I in S1 Text). We did not see as high correlations for “gastroenteritis” search trends as we did with the search trends for the other terms or with HuNoV GII wastewater values (Fig 4B and Fig J in S1 Text), highlighting the importance of terminology. Our results support the use of digital epidemiology as a new tool for disease monitoring, however it has some limitations, including that media can strongly affect search volume (e.g., if an outbreak is reported in the news), and that it has less fine-scale spatial resolution. In particular for wastewater, the spatial resolution of search volume data affected correlations in a similar way as we saw for conventional epidemiological data, where many lag times showed significant correlation. Taken together this suggests that digital epidemiological data on search volume is best used in combination with other data sources such as wastewater and traditional epidemiological surveillance—each source brings distinct advantages and biases, making them particularly useful when considered together.

thumbnail
Fig 4. Correlations of HuNoV GII wastewater values with Norovirus-related search terms.

A. Time course showing daily YC HuNoV GII/PMMoV values (black dots) and normalized search term trends for the norovirus-related terms “norovirus” (orange), “stomach flu” (lavender), and “gastroenteritis” (green) for the Detroit metro region. B. Cross correlations comparing HuNoV GII/PMMoV values from all 5 WWTPs and search term trends for “stomach flu”, “norovirus”, and “gastroenteritis” in the state of Michigan and the Detroit (Det) metro area. Pearson’s correlation coefficients (r) are displayed in a heat map with higher values in blue and lower values in yellow.

https://doi.org/10.1371/journal.pwat.0000198.g004

Conclusions

In this study we assess the comparability of high spatiotemporal resolution HuNoV GII levels in wastewater from 5 WWTPs, to syndromic, outbreak and search trend data over the span of a year. Our results suggest that wastewater monitoring of HuNoV GII leads or concurs with other epidemiological monitoring methods, but correlations between wastewater and other data sources varied by the degree of overlap between the sewershed and the population catchment of the other data source. For example, the cross correlation values obtained when comparing state syndromic data to wastewater HuNoV GII values varied greatly between WWTPs (0.48 to 0.88). The lowest value was seen for the JS WWTP, which collects from a population of only about 90,000 individuals and represents less than 1% of the state population. JS HuNoV GII wastewater values exhibited a sharp peak in early March, unlike many other WWTPs that had elevated levels over a more extended time period. This combination of low population overlap and pronounced outbreak peak in JS likely accounts for the lower correlation values seen for JS. Similar to the state-level syndromic data, NoroSTAT aggregates clinical data across large geographical areas. There is a broad interest in defining the lead time of wastewater data compared to conventional surveillance data. Due to their poor geographic overlap with community wastewater data, the aggregate state syndromic data and national NoroSTAT data are not ideal for assessing the potential lead times of wastewater data. Although the more conventional epidemiological approaches are valuable and can help with forecasting, wastewater-based surveillance can provide a more focused regional picture of the norovirus cases compared to the state data and NoroSTAT that covers such large areas.

Overall, our results suggest that smaller populations and/or closer overlap between the wastewater and syndromic or case populations results in closer temporal correlation. However, there is a limit to how small a sewershed population can be before other factors such as large variability in signals and individual variations in shedding become an issue [40], and this may vary by pathogen.

A limitation of this and future studies is the lack of gold standard case data for norovirus which required us to consider correlations across multiple epidemiological datasets to validate wastewater detection. Additional limitations include that several of our epidemiological data sets could not be resolved at the same spatial geography as our wastewater sewershed data (for example, the Google Trends data could only be resolved to the Detroit region rather than the catchment areas), potentially reducing correlations due to differences in the populations measured, and biasing our lead/lag times if norovirus spatial spread patterns led to one population experiencing increases in transmission patterns before the other population did. It should also be noted that we focused on HuNoV GII in wastewater and HuNoV GI is likely contributing to clinical values in many locations.

Our study adds to the existing literature by considering: multiple treatment plants based in the US spanning a range of population density and urbanicity, high resolution time series data spanning over a year from 600+ wastewater samples, and comparing wastewater data to multiple epidemiological data sets—including taking a digital epidemiology approach [41]. The advantages of wastewater data compared to syndromic and case data include more timely availability and accessibility of data and ability to detect asymptomatic and mild cases who may not seek care. However, these benefits depend on testing frequency, timely reporting using public dashboards, and diversity and representation of the WWTPs being tested. Given that wastewater data for HuNoV GII correlated closely with multiple other syndromic, outbreak, and search term measures of norovirus activity, and did so either concurrently or as a leading indicator, when appropriately implemented wastewater surveillance of HuNoV GII can provide a useful early warning system and an complementary data stream with which to triangulate norovirus patterns—one that does not require healthcare seeking, clinical testing, or inference based on symptom patterns.

Supporting information

S1 Text. Document contains supplementary figures and methods.

Fig A. Comparison of HuNoV GII levels from 5 WWTPs in Michigan from 2021–2022. HuNoV was quantified in influent samples using ddPCR at least weekly and gene copies per liter (gc/L) were plotted over one year. The only exception is TM, where sample collection began later, in January 2022. B and C. Box-and-whisker plots of log10 (B) average weekly gc/L HuNoV values (median, IQR) or (C) HuNoV/PMMoV values for each WWTP for the entire test period. Fig B. Individual graphs of HuNoV GII levels from 5 WWTPs in Michigan from 2021–2022 analyzed by two different methods. HuNoV was quantified in influent samples using ddPCR at least weekly and gene copies per liter (gc/L), shown in orange, as well as HuNoV values normalized to the fecal indicator PMMoV, shown in blue, were plotted over time. The only exception is TM, where sample collection began later, in January 2022. Values from WWTPs in A. Ann Arbor (AA), B. Flint (FL), C. Jackson (JS), D. Tecumseh (TM), and E. area around Ypsilanti Community (YC). Fig C. Seasonal variations in PMMoV-normalized HuNoV GII levels in wastewater. A. HuNoV was quantified in influent samples using ddPCR at least weekly. The only exception is TM, where sample collection began later, in January 2022. The weekly averages of HuNoV/PMMoV were analyzed and a box and whisker plot were used to display the log10 values (median, IQR). B. Seasonal HuNoV GII values (median, IQR) for each individual WWTP are presented using box-and-whisker plots. Note, no fall values for TM were obtained. The definition of the seasons is meteorological, beginning on the 1st day of the equinoxes or solstices. Fig D. Seasonal variations in HuNoV GII levels (gc/L) in wastewater. A. HuNoV was quantified in influent samples using ddPCR at least weekly. The only exception is TM, where sample collection began later, in January 2022. The weekly averages of HuNoV GII GC/L were analyzed and a box and whisker plot were used to display the log10 values (median, ICR). B. Seasonal HuNoV GII values (median, IQR) for all WWTPs are presented using box-and-whisker plots. Note, no fall values for TM were obtained. Fig E. Cross correlations of HuNoV GII wastewater values in gc/L with weekly school reported gastrointestinal illnesses. Graphs show Pearson’s correlation coefficients (r) and probability values (p) determined by comparing weekly average wastewater data in HuNoV II gc/L, to school reported GI illnesses. A. Weekly number of school reported GI illnesses in TM schools was compared to weekly average TM WWTP HuNoV GII values. B. Weekly County-level school reported GI illness values normalized for attendance were compared to weekly average WWTP HuNoV GII values for AA and YC. Cross correlations were tested for a lead (14 or 7 days before), same timing, or a lag (7, 14, 21 and in some cases 28 and 35 days later) in school reported GI illnesses compared to HuNoV GII values. Note: * = p<0.05, ** = p<0.01, and *** = p < .001. Fig F. Correlation of HuNoV GII gc/L wastewater values with weekly GI related emergency department visits at regional and state levels. A. Graph showing Pearson’s correlation coefficients (r) and probability values (p) determined by comparing HuNoV GII wastewater values (in gc/L) from the TM, AA, and YC WWTPs to the weekly total GI related emergency department visits reported for patients from the corresponding counties (Lenawee and Washtenaw), and the TM zip code, for January–April 2022. B. Graph showing Pearson’s correlation coefficients (r) and probability values (p) determined by comparing wastewater data as in A, to the weekly percentage of GI related emergency department visits compared to total visits. Data was reported for patients from the TM WWTP zip code and associated county (Lenawee). State (MI) hospital values were also compared to HuNoV GII wastewater values from all 5 WWTPs for January–April 2022. Cross correlations were tested for a lead (14 or 7 days before), same timing, or a lag (7, 14, and 21 days). Note: * = p<0.05, ** = p<0.01, *** = p<0.001. Fig G. Area contained in google search trends dataset. A. Map showing the lower peninsula of Michigan. B. A zoomed in version of southeast Michigan with the region included in the Detroit Metro google search trend dataset in mauve. Each of our WWTP catchment areas is shown in gray and labeled. Base map tile/data was provided by (c) OpenStreetMap and contributors, CC-BY-SA [openstreetmap.org/copyright]. Fig H. Cross correlations comparing HuNoV GII/PMMoV values from 5 WWTPs and search term trends for “Norovirus” in the state of Michigan (left) and the Detroit Metro area (right). Pearson’s correlation coefficients (r) and probability values (p) were determined. Cross correlations were tested for a lead (-21, -14, -7 days), same timing, or a lag (7, 14, 21, days) in search term trends compared to HuNoV GII/PMMoV values. Note: * = p<0.05, ** = p<0.01, *** = p < .001. Fig I. Cross correlations comparing HuNoV GII/PMMoV values from 5 WWTPs and search term trends for “Stomach flu” in the state of Michigan (left) and the Detroit Metro area (right). Pearson’s correlation coefficients (r) and probability values (p) were determined. Cross correlations were tested for a lead (-21, -14, -7 days), same timing, or a lag (7, 14, 21, days) in search term trends compared to HuNoV GII/PMMoV values. Note: * = p<0.05, ** = p<0.01, *** = p < .001. Fig J. Correlations of WWTP data with “gastroenteritis” search term trends. Correlations were performed to determine the Pearson’s correlation coefficients (r) for HuNoV GII/PMMoV values for all 5 WWTPs (AA, FL, JS, TM, and YC) with search term trends for “Gastroenteritis” (GE) in the Detroit Metro (Det Met) area and Michigan (MI). Table A. Primers, probes, and positive control used in RT-ddPCR.

https://doi.org/10.1371/journal.pwat.0000198.s001

(DOCX)

Acknowledgments

We thank the Ann Arbor WWTP, Flint WWTP, Jackson WWTP, Tecumseh WWTP, and Ypsilanti Community WWTP, Lenawee County Health Department, Jackson County Health Department, Genesee County Health Department, and Washtenaw County Health Department, especially Kristen Schweighoefer, Laura Bauman, and Juan Marquez.

References

  1. 1. Chachu KA, LoBue AD, Strong DW, Baric RS, Virgin HW. Immune mechanisms responsible for vaccination against and clearance of mucosal and lymphatic norovirus infection. PLoS Pathog. 2008;4(12):e1000236. pmid:19079577
  2. 2. Burke RM, Mattison CP, Pindyck T, Dahl RM, Rudd J, Bi D, et al. Burden of Norovirus in the United States, as estimated based on administrative data: Updates for Medically Attended Illness and Mortality, 2001–2015. Clinical Infectious Diseases. 2021;73(1):e1–8. pmid:32291450
  3. 3. Cannon JL, Bonifacio J, Bucardo F, Buesa J, Bruggink L, Chi-Wai Chan M, et al. Global trends in norovirus genotype distribution among children with acute gastroenteritis. Emerg Infect Dis. 2021;27(5). pmid:33900173
  4. 4. Chhabra P, de Graaf M, Parra GI, Chan MCW, Green K, Martella V, et al. Updated classification of norovirus genogroups and genotypes. J Gen Virol. 2019;100(10):1393. pmid:31483239
  5. 5. De Graaf M, Van Beek J, Koopmans MPG. Human norovirus transmission and evolution in a changing world. Vol. 14, Nature Reviews Microbiology. 2016.
  6. 6. Miura F, Matsuyama R, Nishiura H. Estimating the asymptomatic ratio of norovirus infection during foodborne outbreaks with laboratory testing in japan. J Epidemiol. 2018;28(9). pmid:29607886
  7. 7. Shah MP, Hall AJ. Norovirus Illnesses in Children and Adolescents. Vol. 32, Infectious Disease Clinics of North America. 2018. pmid:29406972
  8. 8. Donaldson EF, Lindesmith LC, LoBue AD, Baric RS. Viral shape-shifting: norovirus evasion of the human immune system. Nat Rev Microbiol. 2010;8(3):231–41. pmid:20125087
  9. 9. Yen C, Wikswo ME, Lopman BA, Vinje J, Parashar UD, Hall AJ. Impact of an Emergent Norovirus Variant in 2009 on Norovirus Outbreak Activity in the United States. Clinical Infectious Diseases. 2011;53(6). pmid:21832262
  10. 10. Vega E, Barclay L, Gregoricus N, Williams K, Lee D, Vinjé J. Novel surveillance network for norovirus gastroenteritis outbreaks, United States. Emerg Infect Dis. 2011;17(8).
  11. 11. Yuan K, Huang G, Wang L, Wang T, Liu W, Jiang H, et al. Predicting norovirus in the United States using google trends: Infodemiology study. J Med Internet Res. 2021;23(9). pmid:34586079
  12. 12. Eftim SE, Hong T, Soller J, Boehm A, Warren I, Ichida A, et al. Occurrence of norovirus in raw sewage–a systematic literature review and meta-analysis. Water Res. 2017;111:366–74. pmid:28110140
  13. 13. Kitajima M, Cruz MC, Williams RBH, Wuertz S, Whittle AJ. Microbial abundance and community composition in biofilms on in-pipe sensors in a drinking water distribution system. Science of the Total Environment. 2021;766.
  14. 14. Huang Y, Zhou N, Zhang S, Yi Y, Han Y, Liu M, et al. Norovirus detection in wastewater and its correlation with human gastroenteritis: a systematic review and meta-analysis. Vol. 29, Environmental Science and Pollution Research. 2022. pmid:35048346
  15. 15. Kazama S, Masago Y, Tohma K, Souma N, Imagawa T, Suzuki A, et al. Temporal dynamics of norovirus determined through monitoring of municipal wastewater by pyrosequencing and virological surveillance of gastroenteritis cases. Water Res. 2016;92. pmid:26874777
  16. 16. Rajko-Nenow P, Waters A, Keaveney S, Flannery J, Tuite G, Coughlan S, et al. Norovirus genotypes present in oysters and in effluent from a wastewater treatment plant during the seasonal peak of infections in ireland in 2010. Appl Environ Microbiol. 2013;79(8). pmid:23396337
  17. 17. Wang H, Neyvaldt J, Enache L, Sikora P, Mattsson A, Johansson A, et al. Variations among Viruses in Influent Water and Effluent Water at a Wastewater Plant over One Year as Assessed by Quantitative PCR and Metagenomics. Appl Environ Microbiol. 2020;86(24).
  18. 18. Kazama S, Miura T, Masago Y, Konta Y, Tohma K, Manaka T, et al. Environmental surveillance of norovirus genogroups I and II for sensitive detection of epidemic variants. Appl Environ Microbiol. 2017;83(9). pmid:28213546
  19. 19. Markt R, Stillebacher F, Nägele F, Kammerer A, Peer N, Payr M, et al. Expanding the Pathogen Panel in Wastewater Epidemiology to Influenza and Norovirus. Viruses. 2023;15(2):263. pmid:36851479
  20. 20. Hellmér M, Paxéus N, Magnius L, Enache L, Arnholm B, Johansson A, et al. Detection of pathogenic viruses in sewage provided early warnings of hepatitis A virus and norovirus outbreaks. Appl Environ Microbiol. 2014;80(21). pmid:25172863
  21. 21. Boehm AB, Wolfe MK, White BJ, Hughes B, Duong D, Banaei N, et al. Human norovirus (HuNoV) GII RNA in wastewater solids at 145 United States wastewater treatment plants: comparison to positivity rates of clinical specimens and modeled estimates of HuNoV GII shedders. J Expo Sci Environ Epidemiol. 2023 Aug 7; pmid:37550566
  22. 22. Park HA, Jung H, On J, Park SK, Kang H. Digital epidemiology: Use of digital data collected for non-epidemiological purposes in epidemiological studies. Vol. 24, Healthcare Informatics Research. 2018. pmid:30443413
  23. 23. Bakker KM, Martinez-Bakker ME, Helm B, Stevenson TJ. Digital epidemiology reveals global childhood disease seasonality and the effects of immunization. Proc Natl Acad Sci U S A. 2016;113(24). pmid:27247405
  24. 24. Osuka H, Hall AJ, Wikswo ME, Baker JM, Lopman BA. Temporal relationship between healthcare-associated and nonhealthcare-associated norovirus outbreaks and google trends data in the United States. Infect Control Hosp Epidemiol. 2018;39(3). pmid:29382406
  25. 25. Burnet JB, Cauchie HM, Walczak C, Goeders N, Ogorzaly L. Persistence of endogenous RNA biomarkers of SARS-CoV-2 and PMMoV in raw wastewater: Impact of temperature and implications for wastewater-based epidemiology. Science of the Total Environment. 2023;857. pmid:36240930
  26. 26. Guo Y, Li J, O’Brien J, Sivakumar M, Jiang G. Back-estimation of norovirus infections through wastewater-based epidemiology: A systematic review and parameter sensitivity. Vol. 219, Water Research. 2022. pmid:35598472
  27. 27. Kauppinen A, Miettinen IT. Persistence of norovirus GII genome in drinking water and wastewater at different temperatures. Pathogens. 2017;6(4). pmid:29019921
  28. 28. Flood MT, D’Souza N, Rose JB, Aw TG. Methods Evaluation for Rapid Concentration and Quantification of SARS-CoV-2 in Raw Wastewater Using Droplet Digital and Quantitative RT-PCR. Food Environ Virol. 2021;13(3). pmid:34296387
  29. 29. Kim HY, Kwak IS, Hwang IG, Ko GP. Optimization of methods for detecting norovirus on various fruit. J Virol Methods. 2008;153(2). pmid:18755218
  30. 30. Haramoto E, Kitajima M, Kishida N, Konno Y, Katayama H, Asami M, et al. Occurrence of pepper mild mottle virus in drinking water sources in Japan. Appl Environ Microbiol. 2013;79(23). pmid:24056461
  31. 31. Whale AS, De Spiegelaere W, Trypsteen W, Nour AA, Bae YK, Benes V, et al. The Digital MIQE Guidelines Update: Minimum Information for Publication of Quantitative Digital PCR Experiments for 2020. Clin Chem. 2020;66(8). pmid:32746458
  32. 32. Kageyama T, Kojima S, Shinohara M, Uchida K, Fukushi S, Hoshino FB, et al. Broadly reactive and highly sensitive assay for Norwalk-like viruses based on real-time quantitative reverse transcription-PCR. J Clin Microbiol. 2003;41(4). pmid:12682144
  33. 33. Persson S, Eriksson R, Lowther J, Ellström P, Simonsson M. Comparison between RT droplet digital PCR and RT real-time PCR for quantification of noroviruses in oysters. Int J Food Microbiol. 2018;284. pmid:30005929
  34. 34. Lu X, Wang L, Sakthivel SK, Whitaker B, Murray J, Kamili S, et al. US CDC real-time reverse transcription PCR panel for detection of severe acute respiratory syndrome Coronavirus 2. Emerg Infect Dis. 2020;26(8). pmid:32396505
  35. 35. Wolfe MK, Archana A, Catoe D, Coffman MM, Dorevich S, Graham KE, et al. Scaling of SARS-CoV-2 RNA in settled solids from multiple wastewater treatment plants to compare incidence rates of laboratory-confirmed COVID-19 in their sewersheds. Environ Sci Technol Lett. 2021 May 11;8(5):398–404. pmid:37566351
  36. 36. Decaro N, Elia G, Campolo M, Desario C, Mari V, Radogna A, et al. Detection of bovine coronavirus using a TaqMan-based real-time RT-PCR assay. J Virol Methods. 2008;151(2). pmid:18579223
  37. 37. Zhang T, Breitbart M, Lee WH, Run JQ, Wei CL, Soh SWL, et al. RNA viral community in human feces: Prevalence of plant pathogenic viruses. PLoS Biol. 2006;4(1). pmid:16336043
  38. 38. Arabzadeh R, Grünbacher DM, Insam H, Kreuzinger N, Markt R, Rauch W. Data filtering methods for SARS-CoV-2 wastewater surveillance. Water Science and Technology. 2021;84(6). pmid:34559069
  39. 39. Havumaki J, Eisenberg JNS, Mattison CP, Lopman BA, Ortega-Sanchez IR, Hall AJ, et al. Immunologic and Epidemiologic Drivers of Norovirus Transmission in Daycare and School Outbreaks. Epidemiology. 2021;32(3). pmid:33652444
  40. 40. Fitzgerald SF, Rossi G, Low AS, McAteer SP, O’Keefe B, Findlay D, et al. Site Specific Relationships between COVID-19 Cases and SARS-CoV-2 Viral Load in Wastewater Treatment Plant Influent. Environ Sci Technol. 2021;55(22). pmid:34738785
  41. 41. Tarkoma S, Alghnam S, Howell MD. Fighting pandemics with digital epidemiology. Vol. 26, EClinicalMedicine. 2020. pmid:32864592