National surveillance pilot study unveils a multicenter, clonal outbreak of VIM-2-producing Pseudomonas aeruginosa ST111 in the Netherlands between 2015 and 2017

Verona Integron-encoded Metallo-beta-lactamase (VIM) is the most frequently-encountered carbapenemase in the healthcare-related pathogen Pseudomonas aeruginosa. In the Netherlands, a low-endemic country for antibiotic-resistant bacteria, no national surveillance data on the prevalence of carbapenemase-producing P. aeruginosa (CPPA) was available. Therefore, in 2016, a national surveillance pilot study was initiated to investigate the occurrence, molecular epidemiology, genetic characterization, and resistomes of CPPA among P. aeruginosa isolates submitted by medical microbiology laboratories (MMLs) throughout the country. From 1221 isolates included in the study, 124 (10%) produced carbapenemase (CIM-positive); of these, the majority (95, 77%) were positive for the blaVIM gene using PCR. Sequencing was performed on 112 CIM-positive and 56 CIM-negative isolates (n = 168), and genetic clustering revealed that 75/168 (45%) isolates were highly similar. This genetic cluster, designated Group 1, comprised isolates that belonged to high-risk sequence type ST111/serotype O12, had similar resistomes, and all but two carried the blaVIM-2 allele on an identical class 1 integron. Additionally, Group 1 isolates originated from around the country (i.e. seven provinces) and from multiple MMLs. In conclusion, the Netherlands had experienced a nationwide, inter-institutional, clonal outbreak of VIM-2-producing P. aeruginosa for at least three years, which this pilot study was crucial in identifying. A structured, national surveillance program is strongly advised to monitor the spread of Group 1 CPPA, to identify emerging clones/carbapenemase genes, and to detect transmission in and especially between hospitals in order to control current and future outbreaks.

. Distribution of carbapenemase genes in CIM-positive P. aeruginosa isolates collected within the Netherlands between 2015 and 2017. In the left table, the results from CIM tests and multiplex PCRs for all isolates included in the study are shown. From 1221 P. aeruginosa isolates that were included, 124 (10%) were CIM-positive. In 107 CIM-positive isolates, a carbapenemase gene was also detected. All bla KPC and bla  PCRs were negative, and are therefore not shown. In the right table, NGS results are shown. From the 124 CIM-positive isolates, 112 were sequenced, and the bla VIM-2 allele was discovered in 83 (74%) isolates. VIM Verona Integron-encoded Metallo-beta-lactamase, IMP Imipenem Metallo-beta-lactamase, NDM New Delhi Metallo-beta-lactamase, GES Guiana Extended-Spectrum beta-lactamase, Carba carbapenemase. www.nature.com/scientificreports/ did not yield a PCR product, and were sequenced; three carried a bla GES-5 allele, but a carbapenemase-encoding gene could not be found in the remaining 14 isolates.

Sequencing revealed a large genetic cluster of bla VIM-2 -containing CPPA. NGS was performed on
112 CIM-positive isolates, and on 56 CIM-negative isolates that were matched to CIM-positive isolates based on MML and sampling year. Demographic data provided by MMLs on the patients from which these isolates derived is available in Supplementary Table S1. In some regions, P. aeruginosa was not submitted or found, so isolates from these regions were not sequenced. A complete geographical overview of CIM-positive and CIMnegative isolates that were included is available in Supplementary Figure S1. Genotypic relationships between isolates were determined using wgMLST (Fig. 1). There was a high degree of genotypic diversity, with large allelic distances often exceeding > 3500 alleles between isolates. However, 75/168 (45%) isolates partitioned closely together with a maximum distance of 35 alleles. This group, designated Group 1, comprised isolates that were all CIM-positive (with the exception of one isolate), all belonged to sequence type ST111 and serotype O12, and all but two isolates carried the bla VIM-2 allele; one isolate carried bla  , and in the other, a CIM-negative isolate, no carbapenemase-encoding gene could be identified.
Within Group 1, several genetic clusters could be seen with few allelic differences between isolates. One CIM-negative isolate (also belonging to ST111/O12) was separated from Group 1 by only 95 allelic differences. All other isolates differed from Group 1 by > 3500 allelic differences. Group 1 also comprised isolates from seven Dutch provinces (Fig. 2).
The composition of the integron regions of six Group 1 isolates was reconstructed by combined short-and long-read sequencing. This showed that all four VIM-2-encoding Group 1 isolates contained an identical class 1 integron (Type A) carrying the intI integrase gene, the bla VIM-2 gene flanked by the aminoglycoside resistance www.nature.com/scientificreports/ genes aac(6′)-29a and aac(6′)-29b, an incomplete qacE gene involved in quaternary ammonium compound resistance, and the sulfonamide resistance gene sul1 (Fig. 3). In the CIM-negative Group 1 isolate, the bla VIM-2 and aac(6′)-29b genes were lost from this integron by deletion (Type B). The IMP-13-encoding Group 1 isolate was also shown to carry an integron (Type C) similar to the Type A integron, but contained the bla IMP-13 gene and the aminoglycoside resistance gene aac(6′)-Ib3. Mapping the Illumina reads of the other VIM-2-encoding Group 1 isolates against these reconstructed integrons showed that all carried the Type A integron.
Group 1 isolates were distinct from internationally-derived CPPA isolates. The wgMLST profiles of the 168 sequenced isolates in this study were compared to 260 complete, annotated P. aeruginosa chromosomal sequences from the National Center for Biotechnology Information's GenBank database (Fig. 4). Group 1 is shown in the zoomed in panel; notably, the isolates sequenced in this study (blue circles) were interconnected and not interrupted by a different P. aeruginosa sequence (white circles). However, four isolates with the bla VIM-2 gene were closely related to Group 1, separated by only 3-66 allelic differences. The first strain (accession no. CP016955) was RIVM-EMC4982 from the Erasmus MC University Medical Center Rotterdam, which was the reference isolate used to design the wgMLST scheme for this study. The second strain, Carb01 63 (accession no. CP011317), originated from Maasstad Hospital, another hospital in Rotterdam, the Netherlands. The third strain, PA38182 (accession no. HG530068), originated from a hospital in the United Kingdom, and has been involved in several major outbreaks 14,15 . The fourth strain, PaeAG1 (accession no. CP045739), was a MDR strain aeruginosa isolates that were sequenced in this study, and the number of those isolates that were Group 1, by province. Map was created by importing a screen capture into Adobe Illustrator 2020 and plotting the origin of isolates using the Type-Ned MRSA website (www. typen edmrsa. rivm. nl).

Antibiotic resistance gene profiles and QRDR analysis. ResFinder analyses showed that all 168
P. aeruginosa sequenced isolates carried the beta-lactamase gene bla PAO , the aminoglycoside resistance gene aph(3′)-IIb, and the fosfomycin resistance gene fosA (Supplementary Table S2). All but four isolates carried   www.nature.com/scientificreports/ the chloramphenicol transferase gene catB7. Twenty-six other genes encoding beta-lactamase production, and 27 other genes associated with aminoglycoside resistance, were also found among the 168 isolates. The most prevalent beta-lactamase gene (56%, 94/168) was the bla OXA-50 -like gene bla OXA-395 . This gene was present in all Group 1 isolates, in 29% (11/38) of CIM-positive isolates that did not belong to Group 1, and in 15% (8/55) of CIM-negative isolates that did not belong to Group 1. Both aac(6')-29a and aac(6′)-29b genes were found in all 73 VIM-2-encoding Group 1 isolates; only aac(6')-29a was present in the CIM-negative Group 1 isolate, and neither gene was present in the IMP-13-encoding Group 1 isolate, confirming the integron compositions of Group 1 isolates after read-mapping. Among other isolates, one of either gene was found in three of the 38 CIMpositive isolates that did not belong to Group 1, but neither gene was found in any CIM-negative isolate that did not belong to Group 1. An analysis of quinolone resistance-determining regions (QRDR) for genes gyrA, gyrB, parC, and parE in all sequenced isolates revealed that the gyrA T83I and the parC S87L mutations were found in 100% of ST111 isolates, but were also common among other sequence types; therefore, no unique mutation pattern could be determined for Group 1 isolates (Supplementary Table S3). The combination of T83I and D87N mutations in gyrA were only found in 60% (3/5) of ST175 isolates. The S87W mutation in parC was exclusively found in ST175 isolates.

Discussion
This national surveillance pilot study unveiled that the Netherlands had experienced an ongoing, nationwide, inter-institutional outbreak of a single, clonal genetic cluster of VIM-2-producing P. aeruginosa over a period of at least three years. It was clear from previous reports that several individual hospitals had already recognized this outbreak within their own settings. After the single-center CPPA outbreak reported by van der Bij et al., a surveillance study into 11 hospitals in the Netherlands in 2012 found that the outbreak by CPPA belonging to ST111 was widespread 10 . The surprising results of that study, however, were not followed by the implementation of a structured, national surveillance program. The current study revealed that ST111/O12 has continued to prevail in the Netherlands, and that the outbreak has involved multiple MMLs distributed over a large part of the country.
Genetic clustering of 168 sequenced isolates revealed that almost half (45%, 75/168) showed high genetic similarity, and, with the exception of two isolates, carried the bla VIM-2 allele. Few differences (≤ 35 alleles) were seen between isolates in this genetic cluster, designated Group 1. All VIM-2-encoding Group 1 isolates possessed identical integron structures; coupled with the fact that Group 1 isolates originated from around the country and were submitted by multiple MMLs, we could confirm that Group 1 CPPA have clonally spread throughout the Netherlands. Notably, Group 1 isolates were genetically distinct compared to the other sequenced isolates, exceeding 3500 allelic differences, and to most publicly-available sequences on GenBank. As P. aeruginosa strain Carb01 63 originated from the Netherlands and was isolated in 2012, the authors suspect that this strain also belonged to the outbreak described by this study. Like Carb01 63, Group 1 isolates belonged to sequence type ST111/O12, the predominant P. aeruginosa lineage in Europe 6,7,10 .
P. aeruginosa ST111/O12 clones exhibit high morbidity and mortality in infected patients 6,11,15 . ST111/O12 was first reported in the Netherlands in 2003, more than a decade before the inclusion period of this study, and then again in 2005 9,13 . It is reasonable to suspect that the inter-institutional outbreak described by this study, and the multicenter outbreak described in 2012 by van der Bij et al., are linked, and may have started much earlier than previously anticipated. It is unclear when ST111/O12 was first introduced to the Netherlands, or how country-wide transmission could have occurred. Since it is known that patient referral networks can contribute to the spread of high-risk clones within a country, the transfer of patients between Dutch healthcare institutions most likely played a role in transmission, especially in cases of unnoticed colonization 17 . To date, there have been no studies analyzing the impact of patient transfers between healthcare institutions in the Netherlands. It is highly recommended that patient transfers include accompanying reports on colonization by highly-resistant microorganisms to limit potential inter-institutional transmission. In case of increasing prevalence rates encountered via a national surveillance system, an additional measure could be a national policy to screen patients for CPPA on admission. Screening should especially be performed in patients with a history of recent hospitalization in another healthcare center reporting a CPPA outbreak. Furthermore, CPPA, including ST111/O12 clones, have been shown to reside in the wet niches of hospitals through the formation of biofilm reservoirs 12,14,18 . These reservoirs are persistent, may resist disinfection, and can disperse CPPA to vulnerable patient populations 19,20 , so identifying and limiting environmental sources of ST111/O12 clones in hospitals are of particular interest. ST111/O12 clones have also been found outside of hospitals in wastewater effluents 21,22 ; as the presence of CPPA reservoirs outside of healthcare institutions was not investigated during this study, community-acquired CPPA infections cannot be ruled out.
This study has some limitations. Firstly, no epidemiological data were collected, so no risk factors could be conclusively linked to the outbreak, and transmission events could not be identified. Secondly, MML compliance with submitting isolates was not checked. Thirdly, not all CPPA could be sequenced with the Illumina platform due to budgetary reasons, so it is possible that additional emerging genetic clusters had been missed among the other included isolates; selection bias was limited by sequencing all CPPA isolates received between 2016 and 2017, and sequencing a random selection of CPPA isolates that had been voluntarily submitted in 2015. Finally, long-read sequencing was only performed for six Group 1 isolates.
As a result of this national surveillance pilot study, MMLs have been advised to continue submitting P. aeruginosa isolates to the RIVM for NGS and surveillance. Only isolates with demonstrable carbapenemase production and/or a carbapenemase-encoding gene(s) have been requested. However, structured, national surveillance is still lacking in the Netherlands, since results from current CPPA surveillance are only reported in the Type-Ned database. National surveillance should include alerting the MMLs involved, and supporting epidemiological Scientific Reports | (2021) 11:21015 | https://doi.org/10.1038/s41598-021-00205-w www.nature.com/scientificreports/ investigations into possible transmission routes; this is especially important for new, emerging clones, such as a clone with bla GES-5 . National surveillance should also collect epidemiological data so that risk factors can be assessed. Additionally, the RIVM developed an in-house wgMLST scheme for P. aeruginosa to investigate the clonality of submitted isolates for this study, but recently several other validated schemes were published that may aid surveillance and outbreak investigations [23][24][25] .
In conclusion, the widespread distribution of Group 1 CPPA throughout the Netherlands went unnoticed for a period of at least three years, and this national surveillance pilot study was crucial in identifying the outbreak. Based on previous reports, it is likely that this inter-institutional outbreak started even earlier than previously anticipated. Therefore, the authors strongly recommend the implementation of a structured, national surveillance program in the Netherlands that incorporates wgMLST to monitor the spread of Group 1 CPPA, to identify emerging clones/carbapenemase genes, and to detect transmission in and especially between hospitals to control current and future outbreaks.

Inclusion of isolates.
In 2016, all Dutch MMLs were sent a letter requesting that P. aeruginosa isolates be sent to the RIVM for a national surveillance pilot study. Criteria for submission were that isolates had a minimum inhibitory concentration of > 2 µg/ml for meropenem or > 4 µg/ml for imipenem (as determined by the MML's preferred antimicrobial susceptibility method), and that one isolate per-person-per-year-per-lab was submitted. For each submitted isolate, MMLs were also requested to provide patient age, patient sex, sampling year, sampling site, and MML location in a secured, web-based database called Type-Ned 26 . During an existing surveillance program on carbapenemase-producing Enterobacterales that began before this study, the RIVM had also received P. aeruginosa isolates from several MMLs, so isolates received in 2015 were also considered.
P. aeruginosa isolates submitted between January 1, 2015 to December 31, 2017 were characterized as follows: species level identification was confirmed using MALDI-TOF MS (Bruker Daltonik, Bremen, Germany), carbapenemase production was assessed using the CIM test 27 , and the detection of bla VIM , bla IMP , bla KPC , bla OXA-48 , and bla NDM genes was performed using multiplex PCR with primers and conditions as previously described 26,27 . CIM-positive P. aeruginosa isolates, and a subset of CIM-negative P. aeruginosa isolates matched by sampling year and MML, were subjected to sequencing.

Whole-genome sequence analyses. NGS was performed using Illumina HiSeq 2500 (Illumina, San
Diego, CA, USA), resulting in reads with 125 bases length. De novo assembly was performed using CLC Genomics Workbench v9.5.3 (Qiagen Bioinformatics, Aarhus, Denmark), and contig sequences with a minimum length of 500 bp and at least 30× average read coverage per contig were used for further analyses. QRDR analysis was performed using the sequence extraction tool in BioNumerics v7.6 (Applied Maths, Sint-Martens-Latem, Belgium), in which extracted sequences were translated and aligned to identify coding sequence changes. Sequence types and serotypes were inferred from NGS data in SeqSphere v3.5.0 (Ridom, Münster, Germany) as well as PAst software from the Center for Genomic Epidemiology 28 . Identification of wgMLST alleles was performed in SeqSphere using an in-house wgMLST scheme comprising 6117 core genes and 325 accessory genes based on the fully sequenced and annotated P. aeruginosa strain RIVM-EMC4982 (accession no. CP016955). Allelic distances between isolates were calculated using BioNumerics v7.6. Genes absent in the sequenced isolates were ignored and not counted as allelic differences.
Antibiotic resistance gene profiles were generated by using the ResFinder program v3.2, and the database available from the Center for Genomic Epidemiology website (https:// bitbu cket. org/ genom icepi demio logy/ resfi nder/ src/ master; accessed 08-05-2020) 31 . For resistance gene identification, a 90% identity threshold and a minimum length of 60% were used as criteria.
Ethical statement. The standard administrative procedure for carbapenemase-producing Enterobacterales was used to collect strains and demographic data 26 . Patient identifiers provided by MMLs were encrypted and then stored in the Type-Ned database, ensuring patient privacy in accordance with General Data Protection Regulation. Ethical approval was not required.

Data availability
Sequence data on the isolates in this paper have been deposited in the European Nucleotide Archive under study accession number PRJEB39528 (https:// www. ebi. ac. uk/ ena/ brows er/ view/ PRJEB 39528).