The strains used in this experiment were all isolated from oropharyngeal swabs taken from men who have sex with men (MSM) attending our Sexually Transmitted Infections (STI) clinic in Antwerp, Belgium in 2019 (Laumen et al., 2021b). Nine N. subflava strains containing the msr(D) gene (azithromycin MIC ≥24 mg/L) were used as donor and two N. subflava and one N. gonorrhoeae strains without this gene were used as recipients (MIC <1 mg/L; Table 1). Genomic DNA was extracted using the EpiCentre^® kit. The DNA concentration (ng/μl) was determined using the NanoDrop^® ND-1000 spectrophotometer (Thermo Scientific). One hundred μicroliter of three different donor pools (P1–P3), each containing a mix of three donor DNA extracts of N. subflava (150 ng/μl), were separately mixed with 100 μl (4.0 McFarland) of the mid-log phase growth (6 h) of three recipient strains: (i) N. subflava (ITM_Ns_9/1: azithromycin MIC 3), (ii) N. subflava (ITM_Ns_45/1: azithromycin MIC 6 mg/L and (iii) N. gonorrhoeae (ITM_Ng_38/1: azithromycin MIC 0.19 mg/L; Table 2). Azithromycin concentration of 1.5× MIC was added as a stress factor. Control experiments did not contain azithromycin and/or DNA. The reaction mixtures were plated on blood agar and incubated for 48 h. One colony from each blood agar culture was selected for azithromycin MIC determination E-test gradient strips (bioMerieux, France). All the experiments were conducted at 36°C and 6% CO₂.

The transformation experiment was performed in a NGmorbidostat. The construction, optimalisation and use of the NGmorbidostat have been described in detail elsewhere (Verhoeven et al., 2019; Laumen et al., 2021a). In brief, the NGmorbidostat is a bioreactor that measures bacterial growth via optical density measurements and is used to assess the evolution of antimicrobial resistance (AMR) over time within a constant temperature (35°C–36°C) and CO₂ range (5.5%–6%). In this experiment, we only used the incubator and turbidity measurement functions with the programme MATLAB, to record the growth rate of N. gonorrhoeae (The Math Works, Inc. MATLAB, version R2015b).

The experiment was conducted in four flasks with a total volume of 15 ml in each of gonococcal (GC) broth supplemented with (1%) IsoVitaleX, henceforth referred as GC medium. The conditions were as: (1) 1.5× MIC azithromycin + DNA from msr(D) containing N. subflava, (2) 1.5× MIC azithromycin, (3) DNA from msr(D) containing N. subflava and (4) GC medium (Supplementary Figure 1). To achieve this, firstly we added 200 μl of N. gonorrhoeae (ITM_Ng_21.021 with Azithromycin MIC 1 mg/L; Table 1) at 4.0–5.0 McFarland in all four flasks. After 6 h, the growth curve reached the mid-log phase and 100 μl of HLR-Azithromycin DNA from pool 1 of N. subflava (150 ng/μl) was added to flasks 1 and 3. At the same time point, azithromycin was added to a final concentration of 1.5 mg/L in flasks 1 and 2. After 24 h, 7.5 ml of the old medium was replaced by fresh medium and an additional 100 μl of HLR-Azithromycin DNA from N. subflava (150 ng/μl; flasks 1 and 3) and 1.5 mg/L of azithromycin of the ITM_Ng_21.021 was added (flasks 1 and 2; Table 1). This process was repeated daily for 7 days, after which the azithromycin concentration was increased to 3 mg/L for another 7 days.

Msr(D) was PCR amplified from N. subflava isolate ITM_Ns_41/1 (Table 1) using primers containing a AT-DUS tag to facilitate inter-species transformation to N. gonorrhoeae, forward primer (5’-GAT GCC GTC TGA ACA AAT GAT AAC TGA GG-3’) and reverse primer (5’-GAA TCA ATA CTG ACC AGC GAC-3’). This amplification was carried out as a touchdown PCR: the initial denaturation consisted of 5 min at 95°C, followed by amplification for 10 cycles at 94°C for 30 s, 55°C for 30 s and 72°C for 3 min. The next stage consisted of 35 cycles, lasting 5 more seconds at each cycle, at 94°C for 30 s, 60°C for 30 s and 72°C for 3 min. A final extension step was carried at 72°C for 7 min. The PCR fragment size was analysed on an agarose gel. The concentration of the amplicon was determined using the NanoDrop^® ND-1000 spectrophotometer (Thermo Scientific). The PCR product was used for transformation using the ‘Transformation in plates’ methodology as described above with 100 μl (150 ng/μl) as DNA donor.

Presence or absence of the msr(D) gene in transformant strains were confirmed using quantitative PCR (qPCR). The DNA of the recipient strains was extracted using the EpiCentre^® kit. The primers used to amplify the internal region of the msr(D) (637–934) were as: Forward (5’-GCG GAG GAA AAG CGA AAA C-3’) and Reverse (5’-ACA GAG CCT TAT CCC CAA ATAC-3’). The master mix was composed by 10× EHF PCR buffer (5 μl), 2 Mm dNTPs (7 μl), 5 μM Primer Forward (3 μl), 5 μM Primer Reverse (3 μl), 3.5 U/μl EHF Taq Polymerase (0.5 μl), Rnase free water (21.5 μl) and DNA (10 μl).The qPCR protocol consisted of an initial denaturation stage at 95°C during 5 min followed by amplification for 45 cycles of 94°C for 30 s, 55°C for 30 s and 72°C for 3 min. This step was followed by the final stage consisting of a single cycle of 72°C for 7 min. The specificity of the amplicon was confirmed by conducting melting point analyses.

To test the stability of the transformed N. subflava, a single colony of HLR-azithromycin N. subflava strain (ITM_Ns_27/1; azithromycin MIC of 24 μg/ml) and a single colony of one transformant strain of N. subflava (ITM_Ns_9/1; azithromycin MIC of 256 μg/ml) were subcultured every 24 h in blood agar plates without additional azithromycin for 7 days, similar to the one described in O’Regan et al. (2010). The azithromycin MICs were tested daily on a single colony from each plate with E-tests (Table 3).

The NGmorbidostat was used to compare the growth curves of N. subflava recipient and transformant strain. In a total volume of 15 ml of GC broth supplemented with 1% IsoVitalex (BD BBL™) for each experiment, 100 μl of a 4.0 McFarland suspension in PBS of N. subflava recipient (ITM_Ns_9/1) or N. subflava transformant (ITM_Ns_9/1transformant) strain was added in triplicate. The growth curves were assessed for 18 h, via measurement of optical density every 20 min. Difference in growth curves was assessed via analysis with R (R Core Team, 2019) package ‘growthcurver’ (Sprouffske and Wagner, 2016) with the data obtained from the NGmorbidostat. R was also used to perform the t-test on the samples to confirm or deny or null hypothesis and to obtain the value of p (Supplementary Figure 2).

For whole-genome sequencing (WGS) analysis, the following samples were chosen: (i) DNA recipients after transformation (ITM_Ns_9/1, ITM_Ns_45/1) and (ii) Transformation in morbidostat ITM_Ng_21.021 (Time points 1, 7 and 14 for AZM + DNA and day 14 for the controls). Genomic DNA was extracted using the MasterPure Complete DNA and RNA Purification Kit (Epicentre, Madison, Wisconsin, United States) and suspended in nuclease-free water. Indexed paired-end libraries were prepared using the Nextera XT DNA Library Prep Kit (Illumina, San Diego, CA, United States) and sequenced on an Illumina MiSeq instrument (Illumina, San Diego, CA, United States). Data are available in GenBank: https://www.ncbi.nlm.nih.gov/sra/PRJNA794044. Processed Illumina reads were de novo assembled with Shovill (v1.0.4; https://github.com/tseemann/shovill) which uses SPAdes (v3.14.0) using the following parameters: --trim --depth 150 --opts --isolate (Prjibelski et al., 2020). The quality of the contigs was verified with Quast (v5.0.2; Gurevich et al., 2013) followed by annotation using Prokka (v1.14.6; Seemann, 2014). WGS assemblies of the donor (ITM_Ns_3/2, ITM_Ns_27/1 and ITM_Ns_36/1) and recipient strains (ITM_Ns_9/1, ITM_Ns_45/1 and ITM_Ns_38/1) were available from a previous study by our group and included in the comparative analysis (de Block et al., 2021). BLAST Ring Image Generator (BRIG) was used for genome comparison (Alikhan et al., 2011). Mauve (Darling et al., 2004) was used to align contigs and MEGAX (Kumar et al., 2018) was used to align DNA fragments. Percent sequence identity of DNA fragments was calculated using Muscle (https://www.ebi.ac.uk/Tools/msa/muscle/, version 3.8.31).

Differences were noted in the azithromycin MIC trajectories in the four flasks (Figure 1). The azithromycin MIC of the N. gonorrhoeae recipient increased by day 5 in the flask containing DNA + azithromycin (condition 1) However, the qPCR of msr(D) remained negative in all samples. WGS of samples from day 7 and day 14 revealed a well-known mutation previously linked to macrolide resistance in N. gonorrhoeae: C2611T (Escherichia coli numbering) in the 23S rRNA gene.

WGS of the N. gonorrhoeae recipient in the azithromycin control on day 14 (condition 2) revealed that the recipient acquired the recently described macrolide resistance-associated mutation (RAM) G70D in the 50S ribosomal protein L4 (rplD; Ma et al., 2020; Laumen et al., 2021a).

There was no increase in azithromycin MICs of the N. gonorrhoeae recipient in the DNA control (condition 3) and the growth control (condition 4).

There was no increase in the azithromycin MIC of the N. gonorrhoeae recipient strain after incubation on agar plates for 48 h with a dsDNA fragment containing msr(D) and a DUS. qPCR analysis confirmed that the msr(D) gene was not acquired by the recipient (ct value >30).

Two PCR-confirmed msr(D)-transformant N. subflava strains (ITM_Ns_45/1_P1 and ITM_Ns_9/1_P1; Table 2) were subjected to WGS to identify the exact integration site of msr(D). DNA sequences including the upstream (6,725 bp) and downstream (9,134 bp) region of msr(D) with a maximum total length of 17,803 bp were extracted for further analyses. Donor (ITM_Ns_3/2, ITM_Ns_27/1 and ITM_Ns_36/1), recipient (ITM_Ns_45/1 and ITM_Ns_9/1) and transformant (ITM_Ns_9/1_P1 and ITM_Ns_45/1_P1) DNA sequences were aligned. This alignment revealed the acquisition of the msr(D) gene at the same site (GCATA-acquisition of msr(D)-ATTGA) in the chromosome in both recipients, 32 bp downstream of a DUS sequence (Figure 2). Genome comparison of donor, recipient and transformant revealed that the transformants had acquired a new msr(D)-containing DNA fragment, which originated from the donor, and was not present in the recipient (Figure 3).

A more global alignment conducted in Mauve illustrates the chromosomal organisation around the acquired msr(D) in the transformant ITM_Ns_45/1_P1 compared to the recipient strain (Figure 4). SNP analysis revealed that transformant ITM_Ns_9/1_P1 had taken up a larger DNA fragment than transformant ITM_Ns_45/1_P1 ((Figure 2; Table 4; Supplementary Table 1). The length between the first and last SNP of transformants compared to the recipient strains was 12,033 bp and 5,113 bp for transformants ITM_Ns_9/1_P1 and ITM_Ns_45/1_P1, respectively. The acquired DNA extended from upstream of msr(D) (7,234 bp in ITM_Ns_9/1_P1; 2,883 bp in ITM_Ns_45/1_P1) to downstream of msr(D) (3,335 bp in ITM_Ns_9/1_P1, 766 bp in ITM_Ns_45/1_P1; Supplementary Table 1).

The uptake-fragment of transformant ITM_Ns_9/1_P1 showed high similarity to ITM_Ns_36/1 (97.57% identical to ITM_Ns_3/2, 92.78% to ITM_Ns_27/1 and 99.96% to ITM_Ns_36/1) and the uptake-fragment of transformant ITM_Ns_45/1_P1 showed high similarity to ITM_Ns_3/2 (99.92% identical to ITM_Ns_3/2, 92.41% to donor ITM_Ns_27/1 and 98.43% to donor ITM_Ns_36/1; Table 4). These data suggest that for both transformants, the complete msr(D) containing fragment was taken up from a single (but different) donor in a single transformation event (ITM_Ns_3/2 as the donor for ITM_Ns_45/1 and ITM_Ns_36/1 for ITM_Ns_9/1_P1).