Space-time analysis of COVID-19 cases and SARS-CoV-2 wastewater loading: A geodemographic perspective.

Severe acute respiratory syndrome - coronavirus 2 (SARS-CoV-2) continues to effect communities across the world. One way to combat these effects is to enhance our collective ability to remotely monitor community spread. Monitoring SARS-CoV-2 in wastewater is one approach that enables researchers to estimate the total number of infected people in a region; however, estimates are often made at the sewershed level which may mask the geographic nuance required for targeted interdiction efforts. In this work, we utilize an apportioning method to compare the spatial and temporal trends of daily case count with the temporal pattern of viral load in the wastewater at smaller units of analysis within Austin, TX. We find different lag-times between wastewater loading and case reports. Daily case reports for some locations follow the temporal trend of viral load more closely than others. These findings are then compared to socio-demographic characteristics across the study area.

Introduction

It has been over two years since the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), which causes coronavirus disease 2019 (COVID-19), caused a worldwide lockdown (WHO, 2020). Unlike previous outbreaks of human pathogenic viruses (e.g., SARS-CoV-1 in 2003, H1N1 in 2009, and the Middle East respiratory syndrome [MERS] coronavirus epidemic between 2012 and 2015), the high transmissivity rate of SARS-CoV-2 made this outbreak unique (Fauci et al., 2020; Godri Pollitt et al., 2020; Petersen et al., 2020; Rodpothong and Auewarakul, 2012). At least part of the high transmissivity is attributable to the lag between viral shedding and the onset of physical symptoms (Petersen et al., 2020). This lag, coupled with the substantial hurdles associated with large-scale testing early in the pandemic, made quarantine measures less effective in reducing transmission (Kucharski et al., 2020). Although many countries have implemented (and continue to implement) physical distancing requirements, mask mandates, and restrictions on indoor and outdoor gatherings, more than 532 million COVID-19 cases have been confirmed worldwide (JHU, 2021). However, many scientists believe the reported cases underestimate the actual count due to problems with testing (availability, false negatives), asymptomatic individuals, and a reluctance for testing by some symptomatic individuals (Alwan, 2020; Tanne, 2020; Wu et al., 2020). Alwan (2020) suggests that one of the main factors contributing to testing reluctance is the financial impact of missing work for those without paid time off or sick leave.

The ability to work remotely influences the financial impact of a positive COVID-19 test. For people who cannot work remotely after a positive COVID-19 test, the financial impacts could be devastating, while those who can work remotely likely feel much less of an impact. Remote work is not an option for many industries focusing on service (e.g., restaurants, hotels, janitorial services), construction, or the gig-economy (e.g., Instacart, Uber), among many others. A positive COVID-19 test for such a worker could mean that they go without pay. For most people, losing income is not a viable option. This potential loss means that many individuals choose to work, regardless of COVID-19 status. Compounding these challenges is that employees in these aforementioned industries skew towards lower-income minority groups (Goldman et al., 2021). As a result, recent research finds a disproportionate representation of minority populations in COVID-19 cases, hospitalizations, and deaths (S. J. Kim and Bostwick, 2020; Muñoz-Price et al., 2020; Wadhera et al., 2020).

Due to testing difficulties (e.g., availability), reluctance to test, asymptomatic individuals, and the lag between viral shedding and physical symptoms, it is challenging to accurately identify COVID-19 community case counts (Wu et al., 2020). This uncertainty has made mounting an appropriate response strategy problematic (Noh and Danuser, 2021). To that end, a growing number of researchers are turning to wastewater epidemiology to assess viral loads in a sewershed (Larsen and Wigginton, 2020). This work has primarily focused on developing surveillance tools to test for SARS-CoV-2 RNA (Ahmed et al., 2020; La Rosa et al., 2020), identify ideal sampling points within a sewer network (Balboa et al., 2021; Yeager et al., 2021), and compare SARS-CoV-2 concentrations in wastewater with reported cases (Medema et al., 2020; Randazzo et al., 2020). Results from these studies highlight the value of wastewater in monitoring the prevalence of COVID-19 within a sewershed and the assocaited communities. However, from a response and mitigation perspective, the size of a sewershed makes it difficult to pinpoint specific sub-communities for treatment, testing, and supplies without supplementing sampling at the wastewater treatment plant with sampling of key manholes in the sewershed.

Increasingly, studies are exploring the efficacy of wastewater monitoring to provide actionable insights for COVID-19 response. This includes the recent work by Scott et al. (2021) who utilized gene markers to monitor the spread of COVID-19 on a U.S. college campus, Pillay et al. (2021) who provided a unique perspective on wastewater monitoring in South Africa, and Saththasivam et al. (2021) who were the first to measure SARS-CoV-2 RNA fragments in wastewater in Qatar. This work continues to evolve with some of the most recent studies taking place internationally, such as in Madrid (Lastra et al., 2022). For more information on the state of the art in monitoring we recommend the recent review by Lahrich et al. (2021). In this study, we aim to contribute to this growing body of research by exploring the spatial variation in SARS-CoV-2 loading at wastewater treatment plants and daily reported COVID-19 cases at the ZIP code level in Austin, Texas. To do this, we begin by exploring the spatio-temporal variation of daily COVID-19 cases for each ZIP code falling within two sewershed catchment areas. We then compare the temporal relationship of SARS-CoV-2 wastewater loading (N-gene copies/day) and daily reported COVID-19 cases to uncover the temporal dynamics of viral shedding and case reporting at a more granular spatial scale. Our cross-correlation time-series analysis identifies the presence of lead periods between spikes in SARS-CoV-2 loading and reported cases when considering the aggregate number of cases in the sewershed and wastewater loadings. However, when we disaggregate and apportion cases to the individual ZIP codes based on the population within the sewershed, the lag periods are generally the same, but variation in the strength of correlations between ZIP code case reports and loading levels become clear. The results of this paper support the use of wastewater analysis as a surveillance tool and highlight how case disaggregation and geospatial analysis can enhance surveillance.

Background

The concept of community vulnerability connects to anthropogenically induced or natural hazards – typically large pollution releases, hurricanes, tornados, and other events with the potential to cause significant damage to the built and natural environment. In this context, vulnerability is defined as a community's susceptibility to harm or the potential for loss within the community (Cutter, 1996; Weichselgartner, 2001). Building a deeper understanding of the link between community vulnerability and hazard risk is critical for minimizing the damage caused by disasters – whether anthropogenic or natural.

The ongoing COVID-19 pandemic, by any standard, is a disaster. Millions of people have died (JHU, 2021), the pandemic continues to cripple economies (Bauer et al., 2020), and the ripple effects of the pandemic will manifest for many years (British Academy, 2021; Horowitz et al., 2021). Simply put, the world was woefully underprepared for an event of this magnitude, as evidenced by a general reluctance to act on early information about the virus (Caduff, 2020), the rampant spread of misinformation (H. K. Kim et al., 2020), lack of testing (Caduff, 2020), scarcity of medical supplies (Ranney et al., 2020), and general supply-chain problems (Guan et al., 2020). Moreover, the virus exposed many inadequacies associated with our collective ability to respond to this type of event, highlighting how the most vulnerable communities face the highest risk of COVID-19 illness and death, due, in part, to difficulties associated with monitoring community spread and prevalence (S. J. Kim and Bostwick, 2020).

Socioeconomic health inequities of COVID-19

Post-hoc research on COVID-19 case rates reveals patterns of effects between COVID-19 and socio-demographics that parallel similar research on other environmental and social inequities. The burdens fueled by COVID-19 have fallen disproportionately on vulnerable and marginalized communities. This finding appears to be scale agnostic, with research confirming these patterns across cities (DiMaggio et al., 2020), counties (Liao and De Maio, 2021), states (Karaca-Mandic et al., 2021), and regions (Strully et al., 2021). These burdens range from higher percentages of positive tests to increased hospitalizations and mortality.

More specifically, recent work in Louisiana found that Black Americans made up 70.6% of the COVID-19 hospitalizations and 76.3% of COVID-19 deaths, even though only 36% of the population for the study area is Black (Price-Haywood et al., 2020). Likewise, in Milwaukee, Wisconsin, a cross-sectional analysis revealed that Black males were at higher risk of testing positive for COVID-19 and being hospitalized due to COVID-19 complications than were their racial and ethnic counterparts (Muñoz-Price et al., 2020). Almagro & Orane-Hutchinson (2020) reported similar results, identifying a significant, positive relationship between low-income majority-minority areas and the percentage of positive COVID-19 tests in New York City. Furthermore, Mahajan & Larkins-Pettigrew (2020) and Karmakar et al. (2021) both show the existence of significant racial inequities in COVID-19 mortality using data from across the United States.

Many medical professionals point to the health inequities in the United States as the underlying reason for some of the racial disparities in COVID-19-related health outcomes (Bibbins-Domingo, 2020). For example, pre-existing conditions can make the risk of COVID-19 illness greater. Because minority communities are disproportionately affected by diabetes, hypertension, obesity, and other ailments (Towne et al., 2017; Zarefsky, 2020), their mortality risk for COVID-19 is much higher. In addition to comorbidities, members of minority groups are more likely to be part of the low-wage essential workforce, such as bus drivers, custodians, service workers, and other front-line industry workers (Bibbins-Domingo, 2020). Indeed, Almagro & Orane-Hutchinson (2020) tested and confirmed that working in transportation, construction, and other service industry occupations was positively and significantly associated with positive COVID-19 tests. Because these workers typically do not have paid time off or sick leave (Schneider and Harknett, 2021), they might feel less inclined to be tested or stay home from work following a positive test or symptom onset. As a result, the actual COVID-19 case count might be much higher in such communities. Thus, even though free COVID-19 testing is available, testing is unlikely to be the most reliable means for monitoring COVID-19 prevalence within a community.

Finally, recent research suggests that commercial activity, specifically the movement of people between work and home, is an important explanatory factor in the spread of COVID-19 (Bontempi et al., 2021; Grubesic et al., 2021). Indeed, recent studies have demonstrated an interesting link between COVID-19 transmission in communities where individuals are likely to work at “risky” facilities with the potential for high infection rates (Towers et al., 2022). In short, individuals who commute to and from jobs that place them at higher risk of infection are in-essence two-way vectors for transmission. They may transport the virus to and from the locations where they work and live.

Proactive monitoring

Proactive monitoring is essential for preparedness and effective resource allocation (Lin Moe and Pathranarakul, 2006). Due to the concern regarding underestimation of the actual case count and level of community spread of COVID-19, scholars have been advancing SARS-CoV-2 monitoring in wastewater. Farkas et al. (2020) detailed that viral concentrations in wastewater can provide an essential indicator of viral disease prevalence in a community. While the methods for concentrating and detecting SARS-CoV-2 in wastewater continue to be developed, many researchers have successfully extracted SARS-CoV-2 RNA from untreated sludge samples or other bulk wastewater sampling points (Balboa et al., 2021; Palmer et al., 2021). Quantifying the viral concentration in wastewater samples typically involves reverse-transcription, quantitative, real-time polymerase chain reaction (RT-qPCR), digital PCR, or digital-drop PCR (ddPCR) methods. For example, Ahmed et al. (2020) employed RT-qPCR to determine the SARS-CoV-2 concentration in wastewater to estimate the prevalence of COVID-19 among those who lived within the target Australian sewershed, while (Liu et al., 2020) reported that the ddPCR test returned fewer false negatives in low viral load specimens than did RT-qPCR methods. The results from Suo et al. (2020) also suggest that ddPCR might be better than RT-qPCR for minimizing false negatives.

Wastewater-based epidemiology lends itself to an early-warning and preparedness system for decision-makers and health officials. For example, several recent studies document the presence of SARS-CoV-2 in wastewater weeks before the first case was actually reported (La Rosa et al., 2020; Medema et al., 2020; Randazzo et al., 2020). Moreover, results from other studies suggest that one can use wastewater monitoring to track the prevalence of the virus as it ebbs and flows within a community over time. For example, Ahmed et al. (2020) explored this temporal dynamic and found some initial alignment between the decline in the first COVID-19 case wave and wastewater loading. Kumar et al. (2021) also published promising results, graphically depicting a correlation between wastewater SARS-CoV-2 concentrations and reported COVID-19 cases with a wastewater lead time of 1–2 weeks. Weidhaas et al. (2021) report a similar wastewater lead time.

Interestingly, Nemudryi et al. (2020) documented an even more dynamic relationship that varied depending on a waning (decreasing) or waxing (increasing) of COVID-19 case numbers. During the waning period (decrease in cases after initial surge), patient-reported symptom onset preceded detectable levels of SARS-CoV-2 in the wastewater by eight days,3 with a spike in positive COVID-19 tests two days after wastewater detection. Conversely, during a waxing phase (resurgence), symptom onset was found to precede wastewater detection by five days, with positive tests following four days after SARS-CoV-2 was detected in the wastewater. Given that this analysis was performed in the relatively small city of Bozeman, Montana, the lag between symptom onset and wastewater detection might be longer than expected due to the difficulty in detecting SARS-CoV-2 in wastewater for the small COVID-19-positive population.

These studies highlight the potential for wastewater monitoring as an early-warning system and method for estimating the prevalence of COVID-19 within the communities served by a given sewershed. However, one drawback of current approaches is the lack of explicitly accounting for the number of people or facilities that are contributing to the wastewater. For example, consider the equation to estimate viral gene copies used in this analysis:

Eq. (1) can then be transformed to estimate the number of ill individuals within the sewershed (Weidhaas et al., 2021) with some assumptions about how much effluent each individual contributes. That said, the estimate of total number of infected individuals based on SARS-CoV-2 concentration critically depends on knowing the size of the population contributing wastewater in the sewershed.

Unfortunately, sewersheds rarely align with administrative boundaries; instead, they often bisect locations where the population is known (e.g., census tracts, counties, ZIP codes). This is because many sewer systems are designed for gravity wastewater collection systems (GWCS) or pressure collection. The former relies heavily on the elevation changes across the landscape to control where and how wastewater flows through the system (Islam, 2017). The latter depends on a series of pumps to move the water through the system. To make the system as efficient as possible, they are designed and modified to optimize flow, which will not always conform to administrative boundaries (Miszta-Kruk, 2016). There also are reporting challenges with COVID-19 cases. Agencies typically assign COVID-19 case counts to each administrative unit (e.g., tract or ZIP code). As a result, researchers interested in estimating case counts or positive case rates based on gene copies of SARS-CoV-2 in wastewater must apportion the population and the number of cases to match the service area of the sewershed. If one has information about the population with sufficient granularity, it is possible to geocode (i.e., add latitude and longitude coordinates) household addresses to provide a highly accurate estimate of the population within the sewershed (Weidhaas et al., 2021). However, this level of granularity is not always available, leading to the use of rough population estimates based on large regional areas (Ahmed et al., 2020).

While wastewater-based epidemiology provides an initial indication of COVID-19 community prevalence, it remains challenging to identify and effectively direct intervention resources to sub-communities that might be more vulnerable to COVID-19 and actively experiencing community spread. These challenges are especially acute when the sewersheds bisect administrative units, which is the case for Austin, Texas (Fig. 1 ). This work explores an approach that will help facilitate a more spatially granular assessment of COVID-19 cases and SARS-CoV-2 loading in wastewater. To accomplish this task, we utilize daily reported COVID-19 cases in Austin, Texas, and SARS-CoV-2 wastewater data from the two major sewersheds in Austin (Fig. 1). We begin by calculating space-time clusters of COVID-19 case prevalence across the ZIP codes in Austin, TX. Then, we perform a cross-correlation analysis between SARS-CoV-2 loading in wastewater and daily reported COVID-19 cases for the proportion of each ZIP code that falls within the sewershed of interest. Not only does this allow us to determine the lead time between viral shedding and reported cases, but it also allows us to determine whether the reported cases in some ZIP codes follow the wastewater loadings more closely than others. If so, it would provide a more nuanced understanding of which communities are experiencing more extensive community spread and require more response resources.

Fig. 1

Austin study area with the two wastewater treatment service areas highlighted.

Data and methods

Data collection

We obtained case information from the City of Austin's COVID-19 dashboard (APH, 2021). Reported cases are laboratory-confirmed, official city numbers and are reported at the ZIP code level. We gathered case data from May 2020 through January 2021 for each ZIP code in Austin. We transformed the reported cumulative cases to reflect the new daily cases by subtracting the previous day's cumulative case count (t-1) from the current day's cumulative case count (t).

Wastewater collection

Flow-weighted, 24-hour composite samples of primary clarifier effluent (500–1000 mL) were collected from May 2020 – January 2021 from the two major wastewater treatment plants in Austin, Texas: Walnut Creek (WC) and South Austin Regional (SAR) wastewater treatment plants (Fig. 1). The WC sewershed covers 169.2 sq. miles and processed approximately 53 million gallons per day (MGD) during the study period. The SAR sewershed covers 251.2 sq. miles and processed approximately 43 MGD during the study period (Palmer et al., 2021). Samples generally were collected three times per week at the WC wastewater treatment plant and two times per week at the SAR wastewater treatment plant. Samples were transported on ice to the laboratory at the University of Texas at Austin, mixed, divided into 50-mL aliquots, and stored at −20 °C.

Sample processing and RNA extraction

The wastewater samples were processed and RNA extracted according to the recommended protocol in Palmer et al. (2021). Briefly, triplicate 50-mL aliquots for each wastewater sample were placed in a room temperature water bath for 1.5 h to thaw. Aliquots were then pasteurized at 60 °C for 1.5 h. Two-phase centrifugation at 4500 xg for 5 min at 4 °C resulted in a pellet that was subsequently resuspended in the lysis buffer provided in the MagMax Microbiome Ultra Nucleic Acid Isolation kit (ThermoFisher, Waltham, Massachusetts). Homogenization was conducted in the lysing tubes of the kit via four rounds of disruption in a FastPrep-24 (MP Biomedicals; Santa Ana, California) at 4.0 m/s for 20 s and centrifugation at 13,000 xg for 15–20 s.

RNA was extracted with the MagMax Microbiome Ultra Nucleic Acid Isolation kit according to the manufacturer's instructions utilizing the KingFisher Flex system (ThermoFisher, Waltham, Massachusetts). Extracted RNA was stored at −20 °C. The nucleic acid concentration was measured spectrophotometrically with the Synergy Neo 2 Hybrid Multi-Mode Reader (BioTek; Winooski, Vermont).

SARS-CoV-2 quantification by RT-qPCR

The SARS-CoV-2 concentrations of aliquots were determined via triplicate RT-qPCR reactions on each RNA extract, using a ViiA7 Real-Time PCR System (ThermoFisher; Waltham, Massachusetts) and the CDC nCOV_N2 primer/probe set (Integrated DNA Technologies [IDT], Coralville, Iowa). A standard curve (20,000, 2000, 200, 20, and 2 N gene copies/µL) was prepared via a serial dilution of the 2019-nCoV_N_Positive Control (IDT). RT-qPCR reactions had a total volume of 20 µL, and contained: 5 µL of TaqMan™ Fast Virus 1-Step Master Mix (ThermoFisher; Waltham, Massachusetts), 1.49 µL of the CDC nCOV_N2 primer/probe set mix (for a final concentration of 500 nM of each primer and 125 nM of probe), 8.51 µL of PCR-grade water, and 5 µL of RNA extract. The thermal cycler conditions recommended by the CDC were used: 50 °C for 5 min, 95 °C for 20 s, and 40 cycles of 95 °C for 15 s and 60 °C for 60 s. The limit of detection (LOD) of triplicate N2 assays was determined to be 7.5 N gene copies per reaction via six serial dilutions (100,000 to 1 copy per reaction) of the 2019-nCoV_N_Positive Control SARS-CoV-2 plasmid standard and subsequent 20 replicates of standards with 10, 5, 2.5, and 1 copy per reaction. Because we utilized 50-mL sample aliquots and 200-µL elution buffer, the LOD translates to 6000 N gene copies/L of wastewater. RT-qPCR reactions that did not yield amplification were not assigned a concentration value, and the lack of amplification was noted per the protocol followed by Ahmed et al. (2020a). Samples that did not yield an average concentration of the N gene above the LOD were not included in the correlational analyses. The validity of RT-qPCR results was ensured by including positive controls (2019-nCoV_N_Positive Control, IDT) and negative controls (PCR-grade water) in every RT-qPCR plate. The SARS-CoV-2 concentration in each extract was calculated in units of N gene copies/µL extract by averaging the triplicate reactions for each extract. The average SARS-CoV-2 concentration of each wastewater sample in units of N gene copies/L wastewater was calculated by using the aliquot volume (50 mL), the volume of elution buffer (200 µL), and the concentrations of the extracts. Finally, the average concentration of SARS-CoV-2 is calculated for each wastewater sample and multiplied by the plant's wastewater flowrate for that day (L/d) to obtain the loading of SARS-CoV-2 (N gene copies/d) (Appendix A).

Case data preparation

The area served by each treatment plant and each ZIP code in the Austin area is shown in Fig. 1. It is important to note that several ZIP codes on the periphery of Austin only partially overlap with the sewersheds. Therefore, an essential first step in analyzing the correspondence between SARS-CoV-2 wastewater loading and COVID-19 case count is to apportion cases to each ZIP code based on the percent of the population in each ZIP code that falls within the sewershed area. Several alternative approaches can facilitate this apportionment. First, one could take a purely area-based approach. This approach scales cases by the percent area of each ZIP code that overlaps with the sewershed. However, this simplistic approach fails to account for locations with sparsely (or densely) populated areas. Second, one could use the population of nested administrative units (e.g., block groups or tracts within a ZIP code) that fall within the sewershed to apportion cases. This alternative is more accurate than the area-based approach but will likely face similar issues as the first approach if the sewershed bisects the smaller administrative units, if the administrative unit is sparsely populated, or if the population concentrates in a small area of a larger administrative unit. Third, one could utilize the locations of structures in the administrative unit. Here, the cases are scaled to reflect the proportion of structures within each administrative unit for the sewershed. This approach provides a more accurate representation of population locations while also directly accounting for wastewater sources.

We employed the third approach and began by clipping each ZIP code area to match the spatial bounds of the sewershed. Next, we used a comprehensive structure shapefile to identify the count of structures within each ZIP code and the sewershed boundaries. Next, we divided the number of structures falling within the sewershed portion of the ZIP code ( by the total number of structures in the ZIP code to build a scaling factor for case counts. Specifically, we multiplied this proportion by the total daily case count reported for each ZIP code (2), then, after scaling the cases, we calculated the daily case rate for each ZIP code (or ZIP code portion overlapping the sewershed area) at time for each ZIP code per 10,000 residents (3).

The case rate (3) for each time period () for each ZIP code () was used for the space-time cluster analysis to facilitate a more objective comparison of COVID-19 cases across areas with different populations.

Daily cases – wastewater correlation analysis

Sewershed level

We performed cross-correlations between the wastewater loading at each treatment plant (N gene copies/day) and the apportioned new daily cases4 in aggregate (total for all ZIP codes in the sewershed area) and to each ZIP code individually. Where the wastewater data are concerned, the WC wastewater treatment plant contained three dates (10/13/20, 10/14/20, 10/15/20) with two loading values. For these cases, we took an average of the two values. In addition, loading data from 11/23/20 – 12/2/20 and 12/9/20 – 12/30/20) were omitted from both wastewater datasets because the samples collected on these days had a holding time at 4 °C that exceeded 7 days. We used loading and case data collected between 7/1/20 and 1/11/21 for the SAR sewershed and 5/17/20 to 1/6/21 for WC sewershed.

For the sewershed area analysis, we plotted SARS-CoV-2 loading data for WC and SAR against new daily reported cases. We applied a temporal lag (t + 1, t + 2,…t + 10) to the reported date for the loading data in order to shift it forward in time because SARS-CoV-2 detection in wastewater often precedes symptom onset and case reporting. Next, we used case data corresponding to the date of the lagged loading values and calcualted a Pearson correlation to compare. For example, when analyzing a lag of one day, the COVID-19 case data from 7/2/20 would be paired with the 7/1/20 loading data.

ZIP code analysis

We began the ZIP code level analysis by z-score standardizing both the COVID-19 case counts and SARS-CoV-2 wastewater loadings. We then fit smoothing curves to both datasets using the locally weighted scatterplot smoothing (LOWESS) technique (Nemudryi et al., 2020). The model uses a local weighting value to determine the influence that surrounding points have when determining the value of the smoothed estimate. After several iterations utilizing different weights, we determined that 0.20 offered a good balance between over and underfitting for the wastewater data and 0.166 for the case data. Next, the wastewater data was iteratively lagged relative to the case data in one-day increments beginning at t-15 days and ending and t + 15 days. Next, we calculated a Pearson cross-correlation value between wastewater loading and daily case numbers for each lag period, subsequently recording the lag at which the highest cross-correlation value (max synchrony) was found. Finally, we added both of these measures as an attribute of the ZIP code shapefile for mapping. The lag at which max synchrony occurs provides a general measure of how far ahead wastewater loading preceded case reporting, while the max synchrony value itself will measure how strong the relationship between wastewater loading and reported cases is for each ZIP code.

Space-time analysis

We also performed a space-time cluster analysis of the daily COVID-19 case rates using apportioned case numbers to identify the spatial and temporal clusters of COVID-19 across the sewershed ZIP codes. We employed a space-time variation of the local Moran's I statistic, which determines the percentage of time that each ZIP code was a member of a hotspot (high case numbers surrounded by areas of similarly high case numbers), cold spot (low case numbers surrounded by areas of similarly low case numbers), or outlier (low cases surrounded by high cases, or vice versa) throughout our study period. The local Moran's I is defined as follows:

Where is the attribute for feature , is the global average of attribute , is the weight between feature i and j, and is the variance for attribute . We classify a feature (i.e., ZIP code area) as part of a hotspot, cold spot, outlier, or not significant based on the value of and its associated z-score and p-value. In addition to spatial neighbors , the space-time variation of the local Moran's I considers temporal neighbors ( at a lag of two-time steps, allowing us to evaluate the persistence of clusters across both space and time. We evaluate these results in light of several socioeconomic measures, including a social vulnerability index calculated by the Centers for Disease Control (CDC).5

Results

Once we apportioned the new daily cases to each ZIP code, we conducted cross-correlations between cases and wastewater loading at the aggregate sewershed level and at the ZIP code level. Both of these were performed according to the respective sewershed. When using case data aggregated to the sewershed level, the cross-correlations revealed a significant positive correlation between loading and new daily cases for both sewersheds (WC: r = 0.697, p-value <0.001; SAR: r = 0.648, p-value <0.001). Fig. 2 displays the lag time plot and correlation coefficients. The results indicate a significant positive correlation (p-value < 0.05) for all lag times (Appendix B); however, the correlation coefficient was the highest when wastewater was shifted backwards in time by 1–2-days for WC and 1-day for SAR. Both treatment plants had a second local maxima in correlation coefficient when wastewater data preceded cases by 7- and 8-days (Fig. 2). These maxima align with previous findings, which reported SARS-CoV-2 trends in wastewater appearing 2–8 days before reported COVID-19 cases (Karthikeyan et al., 2021).

Fig. 2

Correlation between SARS-CoV-2 wastewater loading and apportioned aggregate case count for each sewershed at different time lags. Lag time indicates how far wastewater data precedes case data.

We then performed a similar cross-correlation analysis but used the apportioned case data for each ZIP code. Using the ZIP code as the unit of analysis allowed us to map the correlation value for each ZIP code and the lag when max synchrony occurred (Fig. 3 ). The results reveal variation in the strength of the correlation and the lag when the correlation is strongest. Overall, the correlation between wastewater loading and new daily reported cases for each ZIP code was generally high, with an average Pearson value of 0.79 (p < 0.05) and 0.85 (p<0.05) for WC and SAR, respectively.

Fig. 3

The max synchrony value represents the highest correlation between wastewater loading and case count for each ZIP code. The reported time lag indicates the day (lag) associated with the highest correlation between wastewater loading values and case data. Sewersheds are outlined in green and blue.

Both sewersheds share a certain degree of heterogeneity in max synchrony values (right panel Fig. 3). There was a more noticeable east to west trend associated with the max synchrony values in the WC compared to SAR, with the values getting smaller (less correlated) as one moves from the east to the west. In SAR, the synchrony values do not appear to show as strong of an east-west spatial trend. There are also several ZIP codes with low correlation values within Central Austin. This result is not necessarily surprising given the diurnal nature of the population in these areas (i.e., large daytime, commuter populations; smaller nighttime, residential populations) and also the high population of younger individuals (at UT Austin) who may have been asymptomatic and not tested, which may have contributed to an irregular pattern of wastewater loading and case reporting.

The time lag associated with the highest correlation also varied by ZIP code, and interestingly, by sewershed. Recall that we lagged wastewater loading values against case values in 1-day intervals such that positive lag numbers indicate that SARS-CoV-2 wastewater loadings preceded COVID-19 case reporting. ZIP codes in the SAR sewershed generally have a shorter lag period as compared to those for the WC sewershed. Several ZIP codes in the SAR sewershed had a lag of 0; however, there were several other locations in SAR where wastewater loading preceded reported cases by one and four days, and one ZIP code where max synchrony occurred at a 13-day lag (78744). Also of note in the SAR sewershed were two ZIP codes with a negative lag value (78712, 78705). A negative result indicates that COVID-19 cases precede SARS-CoV-2 wastewater loadings. Although this result is possible, the likelihood that it is spurious increases due to these ZIP codes again being associated with the University of Texas at Austin which was a COVID-19 testing hub and has a large young-adult population.

In the WC service area, except for a few areas, most of the ZIP codes show a max synchrony when wastewater loading preceded positive case reports by seven days. The seven-day lag for WC departs from the aggregated analysis (Fig. 2), which found the highest correlation at a lag between zero and two days. That said, it should be noted that the synchrony values for the 1–7 day lag calculated for each ZIP code were separated by roughly a few hundredths in correlation value. It is also worthwhile to recognize the spike in correlation at the seven-day lag point for the aggregate analysis in Fig. 2. This provides some support for the 7-day finding in the disaggregated analysis.

Cluster analysis

After apportioning the daily COVID-19 time-series data, we used case rate (cases/10,000 people) to conduct a space-time cluster analysis. The analysis identified which ZIP codes were associated with COVID-19 hotspots, cold spots, and outliers as well as how long they remained part of each cluster throughout the study period. For illustrative purposes, we focus on hotspots. Fig. 4 displays the percentage of time that each ZIP code was part of a hotspot cluster (if at all).

Fig. 4

Graphical display of the spatio-temporal COVID-19 hotspot clusters for each ZIP code. Coloring corresponds to the amount of time that each ZIP code was a members of a space-time case hotspot denoted as a percent of the total days (320) of reported case data used in the analysis.

In general, the pattern of COVID-19 hotspots follows an east to west trend with values high in the east and low in the west, signifying that ZIP codes in the east were a spatial and temporal hotspot for COVID-19 cases. In addition, many of the ZIP codes identified as part of a COVID-19 hotspot for a large percentage of the study period were the same ZIP codes with large synchrony values illustrated in Fig. 3.

The east-west trend in COVID-19 space-time hotspots is similar to several underlying socio-demographic patterns. To provide more nuance to these patterns, we illustrate the underlying socio-demographic patterns within Austin at the block group level in Fig. 5 , with the ZIP codes outlined in dark gray. The figure highlights the stark east-west trend for most of the city. ZIP codes located in East Austin have a higher level of diversity (racial and ethnic makeup, top left panel), have a higher population of minority groups (top right), and exhibit higher social vulnerability (bottom right) as compared to West Austin. In contrast high median income areas are concentrated in West and central West Austin (bottom left panel). We also assess the spatial distribution of the percentage of the working population that commutes between 25 and 60 min for work (Appendix C) which does appear to follow a similar pattern to that of the synchrony values in Fig. 4.

Fig. 5

Underlying socio-demographic patterns across the Austin, TX metropolitan layer with the sewershed boundaries overlain. The diversity index (top right) provides an indicator of how diverse each ZIP code is with respect to the demographic makeup of the residents. Areas with higher racial and ethnic diversity receive a higher score while areas with more homogeneity with respect to race and ethnicity receive a lower score. Total minority population (top right) measures the size of the minority population within each block group. Median household income for each block group is illustrated in the bottom left, and the social vulnerability index (top right) provides a composite measure of vulnerability based on many variables identified as corresponding to vulnerability to exogenous shocks.

In addition to the visual patterns illustrated in Fig. 5, we tested the correlation between each ZIP code's max synchrony value (the correlation between the SARS-CoV-2 wastewater loading data and new daily COVID-19 cases) and several of the socioeconomic variables and time as a hotspot (Table 1 ). Simultaneous analysis of ZIP codes yields two significant relationships where max synchrony is positively associated with ZIP codes that were spatial and temporal COVID-19 hotspots and the number of commuters. This result is partially expected, given that hotspots of COVID-19 are likely to be major contributors of SARS-CoV-2 levels in the wastewater and that areas with large commuting populations are positively associated with synchrony values, suggesting a relationship between commuting and viral shedding. Median household income, although not significant, shows only a very small correlation for each sewershed, although it is negative for WC (discussed in more detail below).

Table 1

Pearson correlation values between the max synchrony values calculated for each ZIP code from the correlation between wastewater loading and new daily cases and several socioeconomic1 indicators of interest.2

	Time as Hotspot	Diversity Index	Med. HH Inc.	SVI	Commuting population
Combined Max Synchrony	0.3143*	0.2085	0.1629	0.2426	.2828*
SAR Max Synchrony	0.3117	0.0636	0.3513	0.1086	.4109*
WC Max Synchrony	0.5528**	0.5124**	−0.2703	0.5219**	.1783
Significance (p-value): **<0.01, *<0.05

Interestingly, although a visual comparison of the max synchrony values for each ZIP code (Fig. 3) and the socioeconomic variables in the SAR sewershed (Fig. 5) might suggest some significant correlations, statistical testing revealed the contrary; only one of the relationships is significant. Recall that in ZIP codes with greater diversity, the index values are higher. With growing evidence to suggest that minority groups are disproportionately affected of COVID-19, we expected to see a large positive correlation between the max synchrony value and the diversity and minority index. Similarly, even though median household income was not significantly related to max synchrony in the SAR sewershed, the coefficient is positive, which is also unexpected. The coefficient for commuting population is both positive and significant, providing some evidence to suggest that commercial activity, as represented by traveling to and from work, may increase exposure and transmissivity.

The ZIP codes within the WC sewershed broadly exemplify what we expected to see across the region. Except for median household income (which is in the expected direction), all of the coefficients were significant with the expected sign. Specifically, for ZIP codes in the WC sewershed, those with a larger max synchrony value (more association between reported daily case count and wastewater loading) are positively related to the amount of time the ZIP code was a COVID-19 hotspot, a higher diversity index value, and a higher social vulnerability index value.

Discussion

The work presented here demonstrates how the disaggregation and apportioning of COVID-19 cases to smaller units of analysis can reveal more nuanced information regarding community spread when used in conjunction with the level of SARS-CoV-2 wastewater loading measured via RT-qPCR. Specifically, after apportioning cases to each respective ZIP code within the treatment plant sewershed, our analysis revealed that case counts for some ZIP codes follow the peaks and troughs of SARS-CoV-2 wastewater loading more closely than do others. When assessing in which ZIP codes the case counts best match the wastewater loading over time, we identified some spatial patterns that align with more vulnerable communities. We also found that by working at the ZIP code level rather than the aggregate sewershed area, we can forecast new case reports at an earlier time. Clinical testing for SARS-CoV-2 generally happens within 3 to 9 days after symptom onset (Freeman et al., 2021), with viral shedding beginning before symptom onset. Although the aggregated ZIP code data would provide a foreshadowing of case reports 1–2 days in advance, the ZIP code level analysis suggests that some locations wastewater loading values might precede case reports by 7–10 days. There are several implications of these findings worth further discussion.

First, we conducted cross-correlation analyses to examine the effect of disaggregating case data to a smaller unit of analysis (sewershed level vs. ZIP code level). For both analyses, we used the ratio of structures within the sewershed area for each ZIP code to total structures in order to apportion new daily cases to each ZIP code. Both analyses found that SARS-CoV-2 wastewater loading increased before the reporting of COVID-19 cases. However, there were differing temporal lags for the results. Specifically, the highest lag correlation for the aggregate analysis was between 1 and 2 days, while the highest lag correlation for the disaggregate analysis was generally between 3 and 7 days (and earlier for some), depending on ZIP code. There was a noticable second spike in correlation values for the aggregate analysis at the 7 day mark, although not as strong. These results confirm that wastewater loading is a viable approach for surveillance but that the sewershed-level (aggregate data) might not provide sufficient lead time for preparedness. Thus, the use of disaggregate case data could be better for an early-warning system.

Second, many studies have documented the positive relationship between COVID-19 case rate and vulnerable communities (DiMaggio et al., 2020; Strully et al., 2021). We know that both low-income and minority communities have a higher chance of contracting COVID-19 and suffering from COVID-19 related illness due to their higher exposure rates. Unfortunately, it is also true that testing rates are lower in these communities for a variety of reasons. This analysis offers a new perspective on understanding the geography of community spread. Using temporal cross-correlation between SARS-CoV-2 wastewater loading and COVID-19 cases, we identified the ZIP codes with the highest correlation between COVID-19 loading and daily reported cases at the ZIP code level. Our ability to assess the strength of these correlations is related to how we scaled cases. With a deeper understanding of these specific ZIP codes that are contributing more to the wastewater loading, responders can efficiently allocate human resources and response efforts. Moreover, a recent study has identified a link between the probability of contracting COVID-19 from surface water and areas where wastewater treatment plants (max capacity of >10,000 population equivalents) discharge effluent (Wang et al., 2022). The approach taken in this study could further inform this type of surface water analysis by comparing where wastewater treatment plants are discharging effluent and the case rates of the surrounding communities. This would be especially important for understanding whether vulnerable groups were disproportionately exposed to COVID-19 via wastewater treatment discharge.

Third, an important finding from this analysis was the lack of evidence (at the aggregate level) pointing to a correlation between COVID-19 cases in majority-minority communities and wastewater loading. It was only after we investigated the ZIP codes within respective sewersheds that the results matched our expectations. Still, only one of the sewersheds (WC) showed a significant correlation between max synchrony and the socioeconomic variables. These same relationships were not significant for the SAR sewershed or the region more broadly. One explanation for this result concerns the geographic concentration of minority communities in Austin. Comparing the distribution of the minority population and diversity scores at the block group level in Fig. 5 to the distribution of synchrony values in Fig. 4, two important points are revealed. First, in SAR, the concentration of minority groups within any single ZIP code is not as high compared to WC. Second, the max synchrony values associated with each ZIP code in SAR are comparatively dispersed. While there does appear to be a correlation between minority groups in WC, the case data reported at the ZIP code level does not provide us with enough granularity to determine whether the same is true in SAR. This is the likely reason why the SAR correlations in Table 1 are neither significant nor in the expected direction. That said, it is worthwhile to remember that without disaggregating the data and analyzing these relationships at the ZIP code level, we would not have been able to say anything regarding specific relationships between vulnerable communities, COVID-19 loading, and daily reported cases. As a result, the approach of case disaggregation for sewershed analyses still shows promise as a mechanism for more granular hotspot monitoring, preparedness, and mitigation decisions.

Limitations and future directions

One limitation of this work is related to the reporting of COVID-19 cases, where cases are added to the count on the COVID-19 dashboard on the date that test results become available rather than being added to the count on the date that the specimen was collected; thus, the variability with respect to the length of time for analyzing the clinical sample could affect the strength of the correlation analysis performed herein. Where the lag times are concerned, we cannot say for certain why we see some variation in the peak correlation time between wastewater loading and case reports. It may have something to do with the distance that each ZIP code is to the treatment plant but would need more data to be sure. This is a valuable avenue for future research as it may help disentangle when and to what extent specific communities are experiencing an outbreak. Another limitation of this work is the frequency of wastewater testing (two samples per week at SAR and three samples per week at WC). We also recognize that wastewater monitoring is not available nor feasible in all parts of the world, and thus there are some limitations for the generalizability of the approach detailed here. To that extent, researchers should focus on where and how to place sensors within informal sewer systems (those without a treatment plant) in order to gain insight in to where and to what extent COVID-19 is prevalent in these communities

It is also important to recognize the general pattern uncovered in this analysis compared to the location of treatment plants. With the general east-to-west trend of max synchrony values between SARS-CoV-2 wastewater loading and new daily cases and the treatment plants’ location in East Austin, we cannot rule out the possibility that proximity to treatment plants had something to do with the ZIP codes that had the highest correspondence between cases and wastewater loading. Once SARS-CoV-2 enters the sewer system, it begins to degrade (Li et al., 2021); more work is required to understand the effect of signal decay and distance to a treatment plant.

	SAR	WC
Lag	Pearson Coefficient	Pearson Coefficient
0	0.6481**	0.6972**
1	0.8424**	0.6817**
2	0.4129*	0.6988**
3	0.5790**	0.5359**
4	0.5528**	0.3965**
5	0.6966**	0.3513**
6	0.6139**	0.4594**
7	0.7741**	0.5402**
8	0.7636**	0.5049**
9	0.6580**	0.4975**
10	0.5968**	0.4622**
Significance (p-value):*<0.05, **<0.01