Investigation of Antibiotic Resistome in Hospital Wastewater during the COVID-19 Pandemic: Is the Initial Phase of the Pandemic Contributing to Antimicrobial Resistance?

Since the COVID-19 pandemic started, there has been much speculation about how COVID-19 and antimicrobial resistance may be interconnected. In this study, untreated wastewater was sampled from Hospital A designated to treat COVID-19 patients during the first wave of the COVID-19 pandemic alongside Hospital B that did not receive any COVID-19 patients. Metagenomics was used to determine the relative abundance and mobile potential of antibiotic resistant genes (ARGs), prior to determining the correlation of ARGs with time/incidence of COVID-19. Our findings showed that ARGs resistant to macrolides, sulfonamides, and tetracyclines were positively correlated with time in Hospital A but not in Hospital B. Likewise, minor extended spectrum beta-lactamases (ESBLs) and carbapenemases of classes B and D were positively correlated with time, suggesting the selection of rare and/or carbapenem-resistant genes in Hospital A. Non-carbapenemase blaVEB also positively correlated with both time and intI1 and was copresent with other ARGs including carbapenem-resistant genes in 6 metagenome-assembled genomes (MAGs). This study highlighted concerns related to the dissemination of antimicrobial resistance (AMR) during the COVID-19 pandemic that may arise from antibiotic use and untreated hospital wastewater.

Introduction

Antimicrobial resistance (AMR) has become a critical threat to human health globally. Before the novel coronavirus disease (COVID-19) pandemic, bacterial AMR attributed to 1.27 million deaths around the world in 2019.[1] Since the COVID-19 pandemic started, there has been much speculation about how COVID-19 and AMR may be interconnected.[2] Collignon and Beggs[3] argued that the global resistance rates are driven mainly through dissemination and with COVID-19 preventive measures put in place globally, the prevalence of AMR is unlikely to increase except in countries where many other factors, including poor water and sanitation, may see more transmission of antibiotic resistant bacteria (ARB).

Considering that poor sanitation can be a factor that links COVID-19 and threats of AMR, focusing on the prevalence of antibiotic resistant genes (ARGs) and ARB in untreated hospital wastewater may be needed. This is because during the first wave of the COVID-19 pandemic when there were uncertainties and fear of SARS-CoV-2, all patients who tested positive for SARS-CoV-2 were hospitalized for treatment and/or to isolate them from the general community.[4] Antibiotics including fluoroquinolones, cephalosporins, carbapenems, azithromycin, vancomycin, and linezolid were administered to most hospitalized COVID-19 patients.[4] Broad-spectrum antibiotics (i.e., carbapenems, cephalosporins, etc.) were also recommended as a therapeutic treatment against acute bacterial coinfection among COVID-19 patients in Saudi Arabia[5] and the United States.[6] Notably, many studies reported increased antibiotic consumption in hospitals during the pandemic,[7−10] which raised the potential concern of hospitals contributing to subsequent waves of AMR through its wastewater streams as a result of the pandemic. This is especially concerning since there is evidence of a correlation between increased antibiotic usage and ARG abundances in a Scottish hospital’s wastewater.[11] A recent review on the global hospital wastewater treatment scenario further stated that most developing countries often drain the hospital wastewater into municipal wastewater systems or directly discharge them into water bodies without any prior treatment.[12] As poor sanitation remains a dominant factor in driving AMR,[12] there is a need to understand if ARGs and ARB would indeed be increasing in relative abundance within the wastewater streams generated by hospitals during the first wave of the COVID-19 pandemic.

To detect ARGs and ARB, multiple approaches such as cultivation-based methods, quantitative polymerase chain reaction (qPCR), and metaomics can be used.[13] By cultivating for bacterial isolates from clinical specimens, carbapenem-resistant pathogenic Klebsiella pneumoniae and Acinetobacter baumannii were detected among patients admitted into a Parisian and a New Jersey hospital, respectively, during the COVID-19 pandemic.[14,15] In the New Jersey hospital, it was thought that the gain of carbapenem resistance among patients was due to hospital-acquired infection. However, the clinical specimens were only collected at one time point from the patients and were unlikely to provide definitive proof that carbapenem resistance was not acquired prior to their hospitalization. Since the two prior studies that reported the outbreak of carbapenem-resistant bacteria in hospitals, a latter study tracked the prevalence of carbapenem-resistant genes in the wastewater generated by two Finnish hospitals during 9 weeks in Summer 2020.[7] By means of high throughput qPCR, they detected carbapenem-resistant genes blaGES and blaVIM to be prevalent and abundant in the wastewater but did not correlate them to time and incidence of COVID-19 cases. Furthermore, qPCR is a targeted approach and relies on primer sets derived from known gene sequences and hence would not be able to denote temporal variations in ARG sequences that may occur within the hospital wastewater during the COVID-19 pandemic.

In this study, untreated wastewater was sampled from a hospital designated by the Saudi Ministry of Health (MOH) to be one of the only two hospitals (i.e., Hospital A) that provided treatment to COVID-19 patients in Jeddah, Saudi Arabia, during the first wave of the COVID-19 pandemic (April until July 2020). Untreated wastewater was also sampled at the same duration from another hospital (i.e., Hospital B) that did not receive permission from Saudi MOH to receive any COVID-19 patients. Metagenomics was performed to investigate the relative abundance and diversity of ARGs and to determine the correlation of the relative abundance of ARGs with time/incidence of COVID-19. The mobile potential of ARGs by means of their correlation with intI1 (class 1 integron) was further evaluated. Metagenome-assembled genomes (MAGs) associated with ARGs that correlated with time/incidence of COVID-19 were also assessed. The findings from this study would provide insights into what type of ARGs and ARB may be driven by antibiotic usage in a local hospital during the early phase of the COVID-19 pandemic to enter the wastewater streams, hence highlighting specific concerns related to the dissemination of AMR during the COVID-19 pandemic that may arise from untreated hospital wastewater discharge.

Materials and Methods

Hospital Wastewater Sampling

Two hospitals located in Jeddah, Saudi Arabia were monitored from April 22 to July 9, 2020. The first hospital, Hospital A, was designated to provide treatment for only SARS-CoV-2-infected patients during the initial wave of COVID-19, with the number of patients hospitalized per day during the sampling period described previously[16] (Table S1). Untreated wastewater was collected as 1 L grab samples on the morning (approximately 9 AM) of each sampling date from the underground equalization tank. The daily volume of wastewater generated by Hospital A was 750 m3 during the sampling period. The second hospital, Hospital B, serves as a control and is a psychiatric hospital with no COVID-19 patients treated at this hospital and had an average of 75 ± 6 patients during the sampling period. One liter of untreated wastewater was collected as grab samples from the equalization tank of Hospital B on the morning (approximately 9 AM) of each sampling date. The daily volume of wastewater generated by Hospital B was 200 m3 during the sampling period. Grab samples instead of 24 h composite samples were collected during this period as automated samplers were not available on site for both hospitals, and these samplers were not commercially available for immediate deployment during the early phase of the COVID-19 pandemic. While the grab sampling approach only provides a snapshot view of the wastewater collected at that sampling time (9 AM), the limitation of this approach is alleviated in part by sampling over multiple days from April until July 2021. In addition, sampling wastewater in the morning between 8 and 10 AM was previously reported in a separate study to exhibit relatively greater agreement with a composite of 24 h microbial profiles than at other sampling times.[17] All collected samples were stored at 4 °C for not more than 1 week on site at the hospitals before they were transported to the KAUST laboratory and processed immediately as described in section . Daily transferal of samples from hospital to laboratory during the sampling period was not possible as Saudi Arabia imposed a strict movement restriction order at that time. In addition, the storage conditions of 4 °C used on-site were limited by the resources available at the hospital at that time.

Processing of Wastewater Samples for DNA Extraction

250 mL of each sample was individually filtered through a 0.45 μm membrane (Merck Millipore HAWP09000, Cork, Ireland) inside a biosafety cabinet. The membrane with biomass retained on its surface was aseptically cut into strips of about 1 × 3 cm using a sterile scalpel and then placed inside the lysis buffer tube of a PowerLyzer PowerSoil DNA Isolation Kit (MoBio, San Diego, CA). DNA extraction was then carried out with slight modifications to the manufacturer’s protocol. Briefly, lysozyme (1:100, v/v, 100 mg/mL) and achromopeptidase (1:100, v/v, 1 mg/mL) were added to the lysis buffer tube containing the membrane strips, and the tube was then incubated at 37 °C for 1 h before proceeding with DNA extraction based on the manufacturer’s protocol. In total, 53 samples of Hospital A and 21 samples of Hospital B were extracted for DNA. DNA concentration was measured using an Invitrogen Qubit Broad Range assay kit (Thermo Fisher Scientific, Waltham, MA).

Identification of Antibiotic Resistant Genes (ARGs) and Statistical Analysis

All DNA samples were submitted to the KAUST Bioscience Core Lab for shotgun metagenomic sequencing. Paired-end sequencing (2 × 150 bp) was carried out on an Illumina NovaSeq 6000 platform with the same sequencing depth, with an approximately similar sequencing yield of an average 38,000 Mbp per wastewater sample from both hospitals. FastQC (v0.11.2)[18] was used to investigate sequencing quality, and Trimmomatic (v0.39)[19] was used to remove adapter sequences with default settings. Detection of antibiotic resistome from metagenomic short reads was performed by using Bowtie2 (v2.3.1)[20] to map short reads of each sample against CARD (v3.1.0)[21] with MAPQ > 30. The resulting sam files were filtered stringently to include only reads with ≤1 mismatch to the consensus using mapped.py[22] (https://github.com/christophertbrown/bioscripts/blob/master/ctbBio/mapped.py). The fraction of bacterial DNA was obtained by MetaPhlAn3[23] with the metagenome mode. Our calculation of the ratio of bacterial DNA in the data sets showed almost 99% bacterial DNA in all data sets (Figures S1, S2). Thereafter, the general relative abundance of resistome in each sample was evaluated by eq . To address the correlation between the relative abundance of antibiotic resistome and time (and hence the progression of the pandemic), the Spearman’s rank correlation analysis was conducted based on eq .To address the correlation between each ARG with time, the breadth of coverage for each ARG was first determined by eq , and the relative abundance of each ARG (only for those with the breadth of coverage above 80%) was evaluated by eq . To prevent the false positive discovery of beta-lactamase variants, all mapped reads for beta-lactamase (ARO:3000001) were exported and subjected to a secondary functional validation by RGI (v5.2.0) with “perfect and strict hits only” settings. The settings were chosen to ensure the coverage of single nucleotide polymorphism (SNP). This secondary analysis is critical to rule out the false-positive discovery of beta-lactamase variants which often contain single nucleotide mutations.[24] The Spearman’s rank correlation analysis of each ARG with time was conducted based on eq . The p-value for the correlation analysis was adjusted for the false discovery rate by the Benjamini-Hochberg procedure.To investigate the mobility of ARGs, the relative abundance of intI1 (GenBank Accession Number AAQ16665) was evaluated in each metagenomic data set. IntI1 was used to infer mobility because it was suggested as an indicator of horizontal gene transfer and linked to the dissemination of antibiotic resistance in the bacterial community.[25] Briefly, short reads were mapped against intI1 using Bowtie2[20] with MAPQ > 30. The resulting sam files were filtered stringently to include only reads with ≤1 mismatch to the consensus using mapped.py[22] (https://github.com/christophertbrown/bioscripts/blob/master/ctbBio/mapped.py). The relative abundance of intI1 in each metagenomic data set was evaluated based on eq . The Spearman’s rank correlation analysis of intI1 with ARG was conducted based on eq . The p-value for the correlation analysis was adjusted for the false discovery rate by the Benjamini-Hochberg procedure.

Identification of Antibiotic Resistant Bacteria (ARB)

To investigate antibiotic resistant bacteria (ARB) derived from hospitals, metagenomic data sets were applied for genome binning to determine the presence of ARB. Contigs for metagenomic data sets were independently assembled using MEGAHIT (v1.1.4)[26] with default settings.

Long contigs (>1 kbp) were curated for MGE detection. Contigs were first searched for MGEs from a custom MGE database[27] by BLASTN[28] with an e-value of 1e-7 and an identity percentage of 80% (https://github.com/KatariinaParnanen/MobileGeneticElementDatabase). Contigs with MGEs were exported for ARG detection by RGI[21] with “perfect and strict hits only”. Contigs that carried ARGs of interest were curated for visualization using the GGGENES package[29] in R.[30]

Long contigs (>1 kbp) were further binned for genomes on Anvi’o (v7)[31] with a default workflow. CheckM (v1.1.3)[32] evaluated the quality of the metagenome assembled genome (MAG) and filtered out bins with completion > 75% and contamination < 10%. MAGs that passed the above-mentioned criteria were clustered and dereplicated by dRep (v1.4.3)[33] with default settings. Dereplicated MAGs were classified and annotated for their phylogenetic identities by GTDB-Tk (v1.0.2)[34] with default settings. To characterize ARGs within MAGs, genes were detected by RGI (v5.2.0) against the CARD (v3.1.4)[21] database with “perfect and strict hits only” settings. The ARB incidence rate of each sampling day was measured based on eq .The relative abundance of ARB was evaluated to assess the dynamics of ARB over time. Briefly, short reads derived from samples were mapped to MAGs by bowtie2[20] with MAPQ > 30. The resulting sam files were filtered stringently to include only reads with ≤1 mismatch to the consensus using mapped.py.[22] Recovered read mapping files and MAGs were imported into Anvi’o[31] to calculate the ‘mean coverage Q2Q3 per 10 million reads’ (defined as the average depth of coverage across contigs excluding nucleotide positions with coverage values falling outside of the interquartile range and normalized by the number of reads in the metagenome)[35] as the relative abundance of ARBs in samples.

To identify potential horizontal gene transfer (HGT) of ARGs within wastewater bacterial communities, ARG-encoded contigs within ARBs were curated. The protein coding sequence (CDS) coupled with an open reading frame (ORF) and ribosome-binding sites (RBS) on these contigs was detected by prodigal (v2.6.3). ARBs that carried ARGs of interest were curated for visualization using the GGGENES package[29] in R.[30]

Determination of Antibiotic Concentration

To determine the concentration of antibiotics in wastewater of hospitals, seven selected antibiotics, sulfamethoxazole (belongs to sulfonamides); erythromycin (belongs to macrolides), ciprofloxacin and moxifloxacin (belong to fluoroquinolones), penicillin G (belongs to beta-lactam), Meropenem (belongs to carbapenems of beta-lactam) and tetracycline (belongs to tetracycline), were measured by liquid chromatography with tandem mass spectrometry (LC-MS/MS). Details of methodology were included in our Supplementary Methods. To address the correlation between antibiotic concentration and time/ARGs, the Spearman’s rank correlation analysis was conducted based on eq .

Data Sharing

Raw sequence data sets were deposited into the sequence reads archive (SRA) database of the European Nucleotide Archive (ENA) under the study accession number PRJEB49260 with sample accession IDs ERS9122385-ERS9122464.

Results

Relative Abundance of Antibiotic Resistome in Hospital Wastewater Exhibits No Correlation with Time

First, the Spearman’s rank correlation analysis revealed a positive correlation between the sampling time and the number of COVID-19 patients in Hospital A (Spearman’s rank correlation coefficient (SCC) = 0.91, p-value = 3.16 × 10–21). Then, all 53 wastewater samples in Hospital A and 21 wastewater samples in Hospital B were detected with ARGs that matched against the CARD database. Specifically, the relative abundance of antibiotic resistome with respect to the total reads for both hospitals was similar (0.25‰ for Hospital A and 0.18‰ for Hospital B) at the start of the study period (Figure ). However, the average relative abundance of antibiotic resistome in Hospital A was higher with an overall average of 3.79‰ (Figure A) in contrast to 0.16‰ (Figure B) of Hospital B (p-value = 3.47 × 10–23 by Student’s t test). There was, however, no correlation between the relative abundance of antibiotic resistome and time for both Hospital A (SCC = −0.103, p-value = 0.461) and Hospital B (SCC = −0.018, p-value = 0.898). Specifically, the correlation analysis applied to the number of COVID-19 patients in Hospital A with the relative abundance of antibiotic resistome also suggested no correlation (SCC = −0.037, p-value = 0.881).

Figure 1

General relative abundance of resistome in (A) Hospital A and (B) Hospital B. Wastewater samples were collected from April 22 to July 9, 2020. The red line in Hospital A indicates the number of COVID-19 patients in the hospital.

Relative Abundance of Certain ARGs at the Gene-Level Exhibits a Correlation with Time

A gene-level investigation for different ARG categories was further conducted to determine the change in the relative abundance of each ARG with respect to time. A total of 1166 ARGs were detected in wastewater samples of Hospital A, and 50 ARGs (4.29%) had a positive correlation with time (SCC > 0.3, adjusted p-value < 0.05) (Figure S3). The top 10 ARGs with the highest relative abundance were summarized in Table . 189 ARGs (16.2%) identified in wastewater from Hospital A had a negative correlation with time (SCC < −0.3, adjusted p-value < 0.05) (Figure S4), with the top 10 ARGs of the highest relative abundance showing a negative correlation with time, as shown in Table . By comparing results, resistance against tetracyclines, macrolides, sulfonamides, and carbapenems was uniquely observed in positively correlated ARGs and fluoroquinolones as unique drug resistance was observed in negatively correlated ARGs.

Table 1

Top 10 Abundant ARGs with Positive and Negative Correlations with Time in Hospital A

positive correlation with time or number of COVID-19 patients			negative correlation with time or number of COVID-19 patients
ARG name	av relative abundance (‰)	ARG drug class	ARG name	av relative abundance (‰)	ARG drug class
msrE	28.86	pleuromutilins, phenicols, streptogramins, macrolides,lincosamides, tetracyclines, oxazolidinones	aadA25	2.11	aminoglycosides
sul2	9.66	sulfonamides	AAC(6′)-Ib-cr3	2.09	aminoglycosides, fluoroquinolones
VEB-9	1.17	cephalosporins, monobactams	InuB	1.75	lincosamides
AER-1	1.04	penams	QnrB19	1.64	fluoroquinolones
OXA-209	1.01	cephalosporins, penams, carbapenems	bacA	1.63	peptides
AAC(6′)-Ie-APH(2″)-Ia	0.93	aminoglycosides	ArnT	1.62	peptides
OXA-464	0.59	cephalosporins, penams, carbapenems	QnrS2	1.61	fluoroquinolones
tet(39)	0.57	tetracyclines	AAC(3)-IIe	1.56	aminoglycosides
tet(36)	0.27	tetracyclines	TEM-181	1.55	penems, penams, cephamycins, cephalosporins
VIM-4	0.19	penems, penams, cephamycins, cephalosporins, carbapenems	eptB	1.55	peptides

In Hospital B, a total of 74 ARGs were detected, with no ARGs showing a positive or negative correlation with time.

Determining Specific Classes of Beta-Lactamases in Their Correlation with Time

We found a lot of beta-lactamases had a negative or positive correlation with time in both hospitals. Summarizing results of Hospital A, a total of 16 beta-lactamase families were found to be correlated (either positively or negatively) with time (Table ). This includes variants of blaIMP, blaVIM, blaDIM, blaPER, blaVEB, and blaAER which had a positive correlation with time in the wastewater of Hospital A. It was found that these beta-lactamase families belong mainly to 2 groups, either class B beta-lactamases or minor extended-spectrum beta-lactamases (minor ESBL) of class A.[36,37] Notably, class B beta-lactamases (metallo-beta-lactamases) are all carbapenemases,[36] and minor ESBL referred to those ESBLs that were infrequently detected and geographically restricted in contrast to the majority of clinically isolated ESBLs.[37] On the contrary, variants of blaMIR, blaCMH, blaCMY, blaCTX-M, blaSHV, blaTEM, blaKPC, and blaOKP had a negative correlation with time in Hospital A, and they belong to class C beta-lactamases (which were only resistant to cephalosporins)[36] and the majority of clinically isolated ESBLs of class A (Table ).[36−38]

Table 2

Beta-Lactamases (and Their Variants) Detected in Hospital A That Had a Positive or Negative Correlation with Time and/or Number of COVID-19 Patients

beta-lactamases positively correlated with time or number of COVID-19 patients			beta-lactamases negatively correlated with time or number of COVID-19 patients
class B	minor ESBL of class A	carbapenemases of class D	class C	majority of clinically isolated ESBL of class A	non-carbapenemases of class D
IMP (32, 48)	PER (1, 3, 4, 5)	OXA (209, 212, 373, 420, 464, 56, 96, 97)	MIR (14, 16, 2, 23, 3, 6, 9)	CTX-M (101, 103, 108, 11, 114, 117, 139, 142, 144, 15, 155, 156, 157, 163, 164, 28, 3, 33, 52, 54, 55, 66, 69, 71, 79, 80, 82, 88)	OXA (1, 140, 224, 320)
VIM (17, 26, 27, 33, 34, 35, 4, 43)	VEB (9)		CMH (1)	SHV (12, 134, 183, 64, 86)
DIM (1)	AER (1)		CMY (121, 138, 161, 162, 23, 28, 31, 4, 59, 6, 94)	TEM (105, 150, 168, 171, 181, 183, 206, 214, 231, 232, 237, 54, 76, 90)
				KPC (1, 10, 11, 12, 13, 14, 16, 17, 18, 19, 22, 23, 24, 25, 26, 27, 28, 29, 3, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 4, 40, 41, 42, 43, 43, 45, 46, 49, 5, 50, 51, 52, 54, 55, 6, 7, 8, 9)
				OKP (A-11, A-2, A-4, A-8, B-4,B-5, B-6, B-13, B-17, B-19)

In addition, variants of blaOXA had divergent correlations with time or number of COVID-19 patients in Hospital A. Specifically, 12 OXA variants had a positive correlation, but 4 had a negative correlation with time (Table ). The blaOXA belong to class D beta-lactamase, which is a large family with low similarity across its variants.[39] To investigate differences across these OXA variants, all detected OXA variants were clustered based on their amino acid similarity (at the 90% similarity level). Results clustered variants into 4 groups, which are the OXA-48-like group, OXA-1-like group, OXA-211-like group, and OXA-58-like group (Figure ). Additionally, the OXA-48-like group, OXA-211-like group, and OXA-58-like group were carbapenemases[39] that exhibited a positive correlation with time in our study (Table , Figure ). On the contrary, the OXA-1-like group can be inhibited by carbapenem[40] and had a negative correlation with time (Table , Figure ).

Figure 2

Phylogenetic tree for blaOXA variants detected to be correlated with time in Hospital A. Variants were clustered as an OXA subfamily (i.e., OXA-211-like; OXA-58-like; OXA-48-like; OXA-1-like).

Antibiotic Concentration in Wastewater

Antibiotic concentrations in wastewater from Hospital A (Table S2) and Hospital B (Table S3) were determined. Student’s t test revealed significantly higher concentration in Hospital A for sulfamethoxazole (p-value = 0.0104), erythromycin (p-value = 0.00036), ciprofloxacin (p-value = 0.0143), penicillin G (p-value = 0.0111), and Meropenem (p-value = 0.000005) compared to Hospital B. In Hospital A, sulfamethoxazole had a detected average concentration above 0.3 μg/L (Table S2), which was higher than the other antibiotics measured within Hospital A wastewater. There was, however, no correlation between sulfamethoxazole concentration and time. The other antibiotics in Hospital A (conferring fluoroquinolones, beta-lactams, macrolides, tetracyclines) were found in relatively low concentration of less than 0.1 μg/L (Table S2). The Spearman’s rank correlation revealed a positive correlation of Meropenem (i.e., a carbapenem) concentration with time (SCC = 0.71, p-value < 0.05) in Hospital A but not for the other antibiotics in both Hospital A and Hospital B (Tables S2, S3). Additionally, ARGs resistant to the above-mentioned antibiotics were curated, and the correlation between ARGs’ relative abundance and that of antibiotic concentration was determined. There was, however, no correlation found between ARGs and the associated antibiotic concentration in both Hospital A and Hospital B (Tables S2, S3).

Detection of ARGs’ Correlation with intI1

To monitor the transfer potential of ARGs in our samples, the relative abundance of intI1 was investigated in wastewater of both Hospital A and Hospital B. First, the Spearman’s rank correlation analysis suggested no correlation between the relative abundance of intI1 and time (Hospital A: SCC = −0.004 with p-value = 0.348, Hospital B: SCC = 0.075 with p-value = 0.671). However, the relative abundance of intI1 was found to be positively correlated with the relative abundance of antibiotic resistome in both hospitals (Hospital A: SCC = 0.913 with p-value = 0.001, Hospital B: SCC = 0.459 with p-value = 0.002). In the gene-level correlation analysis for Hospital A, 314 ARGs (26.93%) had a positive correlation with intI1 (SCC > 0.3, adjusted p-value < 0.05), and 7 ARGs (0.60%) had a negative correlation with intI1 (SCC < −0.3, adjusted p-value < 0.05) (Figure S5). In the gene-level correlation analysis for Hospital B, no ARGs showed positive or negative correlation with intI1 (Figure S6). Further combining of the gene-level correlation analysis with time and intI1, it was found that 2 ARGs in Hospital A, namely VEB-9 and mphE, had a positive correlation with both intI1 and time (and number of COVID-19 patients) (Figure ). VEB-9 was one of the top 10 abundant ARGs that positively correlated with time in Hospital A (Table ). For Hospital B, no ARGs were found to be positively correlated with both time and intI1, suggesting the low transfer potential for ARGs prevailing in Hospital B wastewater compared with Hospital A.

Figure 3

Spearman’s rank correlation of ARGs with time and intI1 in Hospital A (designated COVID-19 hospital). Each dot represented an ARG detected in Hospital A. The size of the dot represented the overall average relative abundance of each ARG, and the color of the dot represented the class of each ARG. Spearman’s rank correlation coefficient (SCC) > 0.3 was determined as a positive correlation, and SCC < −0.3 represented a negative correlation. Blue and green dashed lines indicated SCC = ±0.3.

Identification of ARGs Associated with MGEs and ARBs

Assembled contigs of metagenomic reads were performed to investigate the co-occurrence of MGEs and ARGs. ARGs of interest (i.e., mphE, VEB-9) that were positively correlated with time and intI1 in Hospital A were determined. We found 9 contigs from 9 wastewater samples that harbored both MGEs and VEB-9 (Figure ). We found 52 contigs from 45 wastewater samples that harbored both MGEs and mphE (Figure S7).

Figure 4

VEB-9 harbored in contigs among the bacterial community of Hospital A. Each arrow represented each CDS (coding sequence). VEB-9 was labeled in red, and mobile genetic elements were labeled in yellow. The intI1 cassette was marked within the dashed box.

Genome binning from the metagenomic reads was performed to investigate the presence of ARB in the hospital wastewater samples. For Hospital A, a total of 1036 MAGs were recovered from 53 wastewater samples. For Hospital B, a total of 1413 MAGs were recovered from 21 wastewater samples. ARGs within MAGs were identified to determine if the MAG is associated with an ARB or not. The incidence rate of ARB on each sampling day was evaluated, and the overall average of ARB incidence in Hospital A (68.8%) was significantly higher (p-value = 6.22 × 10–6) than in Hospital B (47.1%). ARGs of interest (i.e., mphE, VEB-9) that were positively correlated with time and intI1 in Hospital A were determined for their ARB hosts. Through genome-centric analysis, it was found that VEB-9 was carried by 6 MAGs, namely the MAG identified as Tolumomas recovered from the sample on April 24, the MAG identified as Flavobacterium recovered from the sample on June 26, MAGs identified as Prolixibacteraceae and Akkermansiaceae recovered from the sample on June 27, and MAGs identified as Megamonas recovered from the sample on June 28 and subsequently as Fusobacteriales recovered from the sample on July 4. Three MAGs (Tolumonas on April 24, Prolixibacteraceae on June 27, and Fusobacteriales on July 4) were also found on other sampling time points, but they did not encode VEB-9 (Tables S5–S7). We further investigate the distribution of depth (mean coverage Q2Q3 per 10 million reads) of these MAGs across sampling time points, but no correlation was found between depth and time, and the depth of MAGs was not higher when they carried VEB-9 (Tables S5–S7). mphE was, however, not found in any MAGs.

To assess the horizontal gene transfer of VEB-9 across the bacterial community, we further investigate contigs that harbored both MGEs and VEB-9 (Figure ). The intI1 cassette, which included the integron recombination sites (attI and attC),[25,41] was found in contigs from May 27, June 11, and July 4 (Figure S8). This finding revealed VEB-9 was located within the intI1 cassette, reiterating the transfer potential of VEB-9 across the bacterial community through MGEs. In addition, these 6 VEB-9-encoded MAGs were checked for other ARGs, and it was found that all 6 MAGs also carried other ARGs. Specifically, 4 of them carried carbapenem-resistant genes (Table S4), implying coresistance with VEB-9. Coupled with the finding of Figure , it formed a potential transfer flow of VEB-9 across the bacterial community by horizontal gene transfer.

In addition, we focused on the top 10 abundant ARGs that positively correlated with time and number of COVID-19 patients in Hospital A (Table ) to see if they were encoded within any of the ARBs that are also potentially pathogenic. It was found that OXA-209 was contained within a MAG recovered from the sample collected on June 24 and identified to be Shigella flexneri. OXA-464 was contained within 2 MAGs identified to be Arcobacter butzleri and Acetoanaerobium noterae, from sequences recovered from the sample collected on June 20 and July 4, respectively.

Discussion

In this study, Hospital A is a general hospital that has a larger bed capacity and generates a larger volume of wastewater daily than Hospital B. Hospital B is a specialty hospital that focuses on providing treatment for psychiatric patients and drug addicts and hence was not equipped and approved to provide treatment for COVID-19 patients during the early phase of the pandemic. However, an earlier characterization of the wastewater in Hospital B by means of LC-MS/MS revealed the presence of a wide suite of medications including antibiotics within the wastewater.[42] Given the use of antibiotics in both hospitals on a routine basis, the general antibiotic resistome from both Hospital A and Hospital B was therefore at the same relative abundance at the start of the sampling expedition.

However, only Hospital A experienced an increase in the relative abundance of the antibiotic resistome as the pandemic progressed (Figure ). During the first wave of the COVID-19 pandemic when much was unknown about SARS-CoV-2, the use of antibiotics as a therapeutic treatment for infected persons displaying respiratory distress has been to a certain extent heavily influenced by some of the earlier literature touting the benefits of specific types of drugs. For example, the use of azithromycin (a macrolide) in combination with hydroxychloroquine or the use of doxycycline (a tetracycline) was touted to be effective against SARS-CoV-2 during the early stages of the pandemic.[43,44] Darunavir, albeit not an antibiotic but an antiretroviral medication, carries a sulfonamide moiety and was also being considered as a possible treatment for SARS-CoV-2.[45] Although subsequent studies showed that these medications are ineffective against SARS-CoV-2,[46−48] the use of these antibiotics or their related derivatives within Hospital A during the early wave of the pandemic cannot be discounted. There were no published records of specific antibiotic use within hospitals in Saudi Arabia, and both hospitals did not provide usage rates during the sampling period due to confidentiality concerns. However, in many other countries, significant increases in the use of broad spectrum antibiotics (e.g., amoxicillin/clauvanic acid, azithromycin, doxycycline, etc.) were already documented since the pandemic started.[49,50] In fact, a metadata analysis of literature published up to June 2020 indicated that three-quarters of patients with COVID-19 receive antibiotics, with the prescription of antibiotics significantly higher than the estimated prevalence of bacterial coinfection.[51]

Through metagenomics, it was observed that the relative abundance of general antibiotic resistome exhibited no correlation with time in both hospitals. Instead, specific ARGs associated with tetracycline, macrolides, sulfonamides, and carbapenems exhibited a positive correlation with time/number of COVID-19 patients in Hospital A, suggesting that there was a selection of certain ARGs and not the whole resistome per se, that happened during the sampling period. The storage of wastewater at 4 °C for not more than a week on site at Hospital A before the samples were brought into the lab for processing may also have selected for these antibiotic resistant bacteria, possibly because of the selection force incurred due to the presence of their associated antibiotics within the untreated wastewater. However, the concentrations of these tetracyclines, macrolides, and sulfonamides, with the exception of Meropenem (a type of carbapenem), also did not correlate with time, and neither did they correlate with the relative abundance of their associated ARGs. Previous studies that looked at the ARGs in urban sewage also demonstrated a lack of correlation between the abundance of ARGs with the concentration of antibiotics.[52,53] In addition, other studies in livestock production farms also noted that ARGs quantified by means of quantitative PCR increased with an abundance of intl1 and intl2 and continue to persist in high abundances even though the antibiotic concentrations were undetectable within the environmental matrices.[54−56] These collective observations may indicate that although the initial increase in the relative abundance of ARGs as seen in Hospital A was due to higher antibiotic use than in Hospital B, the subsequent increase in the relative abundance of specific ARGs may be due to a combination of sustained antibiotic use and/or horizontal gene transfer among the bacterial community. Alternatively, as both hospitals have sewer systems that do not separate gray water (e.g., wastes generated from kitchen, lobby, cleaning water, etc.) from black water (i.e., wastes generated from fecal discharge), it is also possible that the increased cleaning frequency and other SARS-CoV-2 sterilization strategies put in place at Hospital A may be contributing to the coselection of selected resistant genes in its wastewater.

However, the role of horizontal gene transfer in disseminating ARG can be further observed in this study since there was a positive correlation among the relative abundance of ARGs, intI1, and time during the first wave of the COVID-19 pandemic. Using beta-lactam resistant genes as an example, Vietnamese extended-spectrum beta-lactamase VEB-9 was identified in Hospital A wastewater to exhibit a positive correlation with intI1 and time. No ARGs were found to have a positive correlation with both intI1 and time in Hospital B wastewater, reiterating that the therapeutic regimen provided to COVID-19 patients in Hospital A may have selected for mobile ARGs in the wastewater. VEB and its variants have been described among various Gram-negative pathogens of multiple genera such as Vibrio spp., Pseudomonas aeruginosa, Acinetobacter baumannii, etc.[57−60] In earlier instances where they were identified, they were usually present in gene cassettes within class 1 integron structures.[61] In this instance, the VEB-9 genes were also found in contigs that contain mobile genetic elements such as transposase and integrase (Figure ). The contigs were associated with MAGs that were identified to include both Gram-negative (i.e., Akkermansiaceae) or Gram-positive bacteria (i.e., Megamonas and Fusobacteriales), indicating that VEB-9 may possibly have broad host ranges during the short duration of April until June 2020. From our current metagenomics approach, none of the MAGs identified with VEB, and the associated mobile genetic elements were of clinical relevance as they are not known to be pathogens and do not contain any virulence factors (data not shown). However, all of the 6 MAGs identified with VEB were present with other ARGs, including 4 MAGs that had carbapenem-resistant genes, suggesting difficulty in relying on different generations of beta-lactams to eradicate these bacterial groups.

Furthermore, MAGs associated with Shigella flexneri and Arcobacter butzleri, both of which are pathogens that can cause gastrointestinal distress upon infection, were detected.[62,63] In particular, Shigellosis is a diarrheal disease that can cause severe manifestations among at risk groups.[64,65] Both MAGs contained carbapenem-resistant genes (OXA-209 and OXA-464) and were recovered from sequences of wastewater samples collected on June 20 and July 4, respectively. These two sampling dates were toward the late stage of our sampling period where we noted that both OXA-209 and OXA-464 genes (classified as carbapenemases of Class D beta-lactamase) and carbapenem-resistant Class B beta-lactamases exhibited a positive correlation with time. In contrast, the majority of the beta-lactamases that exhibited negative trends is resistant to cephalosporins or are ESBLs that remain susceptible to carbapenem (Table ). This observation is in agreement with that of earlier studies from Finland, France, and the US that noted the outbreak of carbapenem-resistant bacteria and the prevalence of carbapenem-resistant genes in hospitals during the COVID-19 pandemic.[14,15,66] The increase in carbapenemase within the hospitals cannot be easily explained as many factors, for example, the potential increase in the use of last-resort antibiotics over time, lapses in core infection prevention, control practices as number of patients surged, etc., could possibly contribute to its increase.

As this study only focuses on tracking the changes in ARG and ARB diversity and the relative abundance within the hospital wastewater during COVID-19, a main limitation is that we were unable to assess if the key findings noted here would apply to that experienced in an outpatient environment. In countries with a strong antibiotic stewardship program, antibiotic use in the community has declined during the pandemic.[67,68] However, in other countries, commonly available beta-lactams like amoxicillin/clauvanic acid remain easily accessible and can potentially be utilized by the general community to evade hospitalization and strict quarantine rules, therefore potentially accounting for the increasing use of this antibiotic.[69] Likewise, retail- and hospital-based pharmacy purchases of antibiotics during April-August 2020 increased relative to prepandemic April-August 2019.[70] Our findings suggest that resistant genes associated with these classes of frontline antibiotics (e.g., tetracycline, macrolide) did increase during the first wave of the COVID-19 pandemic but may not be picked up by the current surveillance systems of AMR (e.g., GLASS[71]) since these systems focus mainly on pathogens and antimicrobial agents such as ESBL and carbapenems.[71,72] Furthermore, those that showed a negative correlation with the COVID-19 pandemic in our study are the ARGs associated with the ESBLs that are geographically widespread and dominant, for example, CTX-M, SHV, TEM, and KPC, which would also have been the focus of the current surveillance systems of AMR. Hence, a targeted approach by means of cultivation and/or qPCR would have resulted in the interpretation that the abundance of those related genes decreased with time during the first wave of the COVID-19 pandemic. Future studies should therefore consider mapping prevalence and abundance of ARGs and ARBs in a nontargeted manner using metagenomics within the urban sewage to determine the impact of the COVID-19 pandemic on community AMR.

In summary, this study showed that during the first wave of the COVID-19 pandemic in Jeddah, Saudi Arabia, there was a positive selection of certain ARGs within the hospital wastewater, and the selection was likely due to the prevailing therapeutic treatment given to patients. Over time and/or with an increase in the number of COVID-19 patients, there was a further positive selection of ARGs resistant toward carbapenem, including those that are associated with mobile genetic elements. This can see potential dissemination of ARGs resistant to last-resort antibiotics by means of HGT. Our metagenomics approach further denoted that the presence of these carbapenem-resistant genes can be associated with contigs of a broad range of pathogenic and nonpathogenic bacteria. It is therefore important to consider enforcing stewardship of proper antibiotic use within both hospitals and the community during the pandemic and also to ensure efficient treatment of hospital and/or municipal wastewater during the COVID-19 pandemic. Collectively, this two-pronged approach can serve to mitigate an unintentional contribution toward threats of AMR during the COVID-19 pandemic.