Enrichment of antibiotic resistance genes within bacteriophage

Introduction The increased use of antimicrobials and development of resistance is at an alarming rate (Aslam et al., 2018). The emergence and spread of antimicrobial resistance has been listed by the WHO as one of the major threats to global health security (Roberts and Zembower, 2021). Specifically developing countries have high prevalence of multidrug resistant bacteria (Ayukekbong et al., 2017). Tanzania, being a developing country, has rapidly grown its economy during the last decade (Mrema et al., 2015) but has also seen a significant rise of antimicrobial resistance in particular due to the inappropriate use of antimicrobials in both the human and animal sectors (WHO Joint External Evaluation of IHR Core Capacities of the United Republic of Tanzania, 2017). Recent prevalence studies of antibiotic usage in Tanzania documented a high usage for children and patients admitted to surgical and pediatric wards (Seni et al., 2020). Such prescription was not based on culture-testing or antimicrobial susceptibility testing, but was commonly prescribed without any tests. Similarly, an excessive usage of third-generation cephalosporins (ceftriaxone) among hospitalized patients could be detected in 51.1% (322/630) individuals (Seni et al., 2020). The emergence and spread of multidrug resistant bacteria to both first- and second-line drugs (e.g., amoxicillin, chloramphenicol, trimethoprim-sulfamethoxazole, extended-spectrum cephalosporins, and fluoroquinolones) is well documented in Tanzania (Moremi et al., 2016). The rate of extended spectrum beta-lactamase (ESBL) producing Escherichia coli and MRSA has been increasing at an alarming rate (Moremi et al., 2014; Katale et al., 2020), with high prevalence of also other resistance genes [blaCTX-Ms (45.7%), SCCmec type III (27.3%), IMP types (23.8%); Katale et al., 2020], with resistance to third generations cephalosporins ranging from 26 to 100% in parts of Tanzania (Masinde et al., 2009; Mushi et al., 2014; Ampaire et al., 2016; Moremi et al., 2016; Sonda et al., 2019). Likely these numbers are however too low due to lack of appropriate microbiological diagnostic capacity and infrastructure in most clinical settings, as well as inadequate resources to implement national wide surveillance programs (Katale et al., 2020). Usage of antibiotics is however not only adding selective pressure for resistant bacteria but has also been shown to be a key factor for induction of temperate bacteriophages and spread of resistance through transduction (Stanczak-Mrozek et al., 2017). An earlier in vivo study in mice by Modi et al. demonstrated the impact of antibiotics at subclinical concentrations to induce bacteriophages carrying resistance genes, able to transduce sensitive bacteria (Modi et al., 2013). Similar effects have been detected in humans undergoing antibiotic treatments with increased abundance of bacteriophages carrying resistance genes (Abeles et al., 2015). However, the latter study was very limited in number of patients included (n = 4). There is therefore a need to investigate such putative induction in higher resolution. We recently described how the presence of specific protozoa (Entamoeba gingivalis) can affect the oral microbial diversity (Stensvold et al., 2021). Herein, we study alterations of microbial diversity and prevalence of antibiotic resistance carrying bacteriophages in a Tanzanian patient group undergoing antibiotic treatment. Materials and methods A cross-sectional study was carried out in June 2019 at Tanga Regional Referral Hospital, Tanzania. The study was approved by the Medical Research Coordinating Committee of the National Institute for Medical Research (NIMR MRCC; reference number, NIMR/HQ/R8.a/Vol.IX/3079). Inclusion criteria for individuals entailed age (18 years or above), no usage of antibiotics during the last 3 months, a medical condition treated with oral antibiotics for at least 3 days (patients only) and with full mental capacities. Individuals were excluded from the study if they had non-infectious immune-modulating sicknesses. All individuals signed a written informed consent before enrolment in the study. The study enrolled 25 patients, who had been prescribed oral antibiotics (3–10 days) for treatment of non-oral infections, and 26 patients at the hospital who had not received antibiotics (e.g., non-infectious reasons for hospitalization for 3–10 days). Saliva samples (5 ml) were collected on day 0 (hospitalization day) and day 3. For all samples, information on age, gender, and treatment was available. Bacterial and bacteriophage DNA isolation Bacterial DNA was isolated using Norgen’s Saliva DNA preservation and isolation kit (Norgen Biotek Corp., Thorold, ON, Canada). For bacteriophage DNA preparations, saliva samples were centrifuged (10,000 g 10 min), and the supernatant passed through sterile filters (0.22 μm) to remove bacterial contaminants. The samples were then processed using a phage DNA isolation kit (Norgen Biotek), according to manufacturer’s instructions, including addition of DNAse I and heat-inactivation thereof. DNA was quantified with Qubit fluorimeter (Life Technologies, Carlsbad, CA, United States). 16S/18S amplicon-based microbiome analysis Amplicon-based microbiome analysis was done as previous described (Ring et al., 2017). Library preparation was performed by Nextera XT DNA Library Preparation (Illumina inc., San Diego, California, United States), and Illumina sequencing was performed on the Hiseq system (Illumina) according to the manufacturer’s instructions. DNA was amplified using a two-step PCR using custom 341F/806R primers targeting the V3-V4 16S regions, and three primer sets targeting the hyper-variable regions V3-V4 of the 18S rDNA gene, and amplicons were sequenced on the Illumina MiSeq (Illumina) using the V2 Reagent Kit. Detection of resistance genes by metagenomic long-read sequencing DNA was prepared for sequencing using Oxford Nanopore Technologies’ Rapid PCR Barcoding Kit (SQK-RPB004) with the following modifications to the manufacturer’s instructions: Double volume of template DNA and “FRM,” as well as 25 PCR cycles instead of 14. DNA libraries were sequenced in R9.4.1 flow cells (FLO-MIN106) in a MinION (Oxford Nanopore Technology) connected to a MinIT with MinIT Release 19.12.5 (MinKNOW Core 3.6.5, Bream 4.3.16, and Guppy 3.2.10). Raw reads were basecalled with the “Fast” configuration of the algorithm. Basecalled reads were analyzed with the mapping tool KMA (K-Mer Aligner; Dillip et al., 2018) version KMA-1.2.22 with the following parameters adapted for nanopore data: “-mem_mode -mp 20 -mrs 0.0 -bcNano -and.” KMA were used for mapping against the following databases https://www.arb-silva.de/no_cache/download/archive/release_132/Exports/Archaea (16S & 18S rRNA from bacteria, archaea, and eukarya), https://www.cbs.dtu.dk/public/CGE/databases/KmerFinder/version/20190108_stable/ (Fungi, plasmids, and protozoa), https://www.cbs.dtu.dk/public/CGE/databases/KVIT/version/20190513/ (viruses), as well as ResFinder database (2020-04-08) and PlasmidFinder (2020-04-02). Dehumanization of the raw reads was performed using Minimap2 against the hg38 human genome reference. Reads unmapped to the human genome were extracted using samtools 1.9 with parameter -f4. All data were uploaded to the ENA browser under accession PRJEB55897 (primary) and ERP140841 (secondary), and can be found in Supplementary File 1. Bioinformatics analysis of sequence data Bioinformatics was done using BION, a newly developed analytical semi-commercial open-source package for 16S rRNA and other reference gene analysis, classifying mostly to species level. The pipeline accepts raw sequence and includes steps for de-multiplexing, primer-extraction, sampling, sequence- and quality-based trimming and filtering, de-replication, clustering, chimera-checking, reference data similarities, and taxonomic mapping and formatting. Non-overlapping paired reads are allowed for analysis, and BION is often accurate to the species level. Statistics of sequence data Analysis of microbiome composition was performed in R version 4.0.32020-10-10 using the packages phyloseq v.1.24.2 and vegan v.2.5-2. Figures were created using ggplot2 v.3.2.0 and plotly v.4.8.0. In the bar plot, taxa were merged to genus level by agglomerating counts within each genus. Alpha-diversity of samples as well as relative abundances of individual genera were compared between groups with Mann–Whitney rank sum tests and adjusted for multiple testing using Bonferroni correction. Differences between groups were assessed with bar plots and Principal Coordinates Analysis (PCoA)
Read More