P. yeei strains CCUG 13493, CCUG 17731, CCUG 32052, CCUG 32054, CCUG 46822, and CCUG 54214 were purchased from the Culture Collection of the University of Gothenburg (CCUG) (Sweden). Strain LM20 was isolated from black shale collected within the Lubin underground copper mine in Poland (Dziewit et al., 2015). P. yeei CCUG 13493R, CCUG 17731R, CCUG 32052R, CCUG 32054R, CCUG 46822R, CCUG 54214R, and LM20R, rifampicin-resistant derivatives of the respective wild-type strains, and Escherichia coli TG1 (Sambrook and Russell, 2001), were used as plasmid recipients. E. coli DH5α (Hanahan, 1983) was the host strain of helper plasmid pRK2013, used in triparental matings (Ditta et al., 1980) (Supplementary Table S1).

All strains were cultured in lysogeny broth (LB) medium (Sambrook and Russell, 2001) at 37°C for E. coli or 30°C for the other strains. Liquid cultures were grown under shaking conditions. When required, the medium was supplemented with 10% sucrose and the following concentrations of selective antibiotics: kanamycin, 50 μg/ml for all strains except Achromobacter sp. LM16R (300 μg/ml) (Dziewit et al., 2015); rifampicin, 50 μg/ml; tetracycline, 2 μg/ml for P. yeei and 20 μg/ml for E. coli strains.

The response to various growth conditions (temperature, pH and salinity) of P. yeei strains was determined as described previously by Dziewit et al. (2013). Motility was tested according to the method previously outlined by Dziewit et al. (2015). The hemolytic ability of different strains was tested by growth on blood agar plates that were incubated at 30 and 37°C under aerobic conditions—the results were read after 24 and 48 h of incubation. Siderophore production was assessed using the chrome azurol S (CAS) agar plate method, as described previously (Dziewit et al., 2015). The plates were incubated in the dark at 30°C for 72 h, after which halo formation around the colonies was assessed. Each isolate was phenotypically characterized using (i) the API 50 CH test system (to differentiate strains based on their ability to metabolize different carbohydrates) and (ii) the API ZYM test system (to detect and identify enzymatic activities); according to the recommendations of the supplier (bioMerieux, Marcy l'Etoile, France). Ureolytic activity was detected by cultivation on Christensen's urea medium, containing urea and a pH indicator (phenol red) (Christensen, 1946) or using rapid urease tests. Proteus vulgaris and/or Klebsiella pneumoniae strains (University of Warsaw collection) were used as positive controls, and P. aminophilus JCM 7686 and/or P. aminovorans JCM 7685 as negative controls (Urakami et al., 1990). The minimal inhibitory concentrations (MIC) of antimicrobial agents were determined for P. yeei strains by using Mueller Hinton II agar diffusion with E-tests as recommended by the supplier (bioMérieux). E. coli ATCC 25922 was used as a control. To confirm reproducibility, antimicrobial susceptibility testing was repeated at least twice. Since no interpretative guidelines for P. yeei susceptibility have been published, we used the criteria proposed by Arias and Clark (2019). MIC values for selected heavy metal ions were established on titration plates using a procedure described previously (Dziewit et al., 2013). The following heavy metal salts were used to prepare appropriate stock solutions in water: NaAsO₂; 3CdSO₄ × 8H₂O; CoSO₄ × 7H₂O; K₂Cr₂O₇; CuSO₄; HgCl₂; NiCl₂ × 6H₂O; NaO₃V; ZnSO₄ × 7H₂O. Each microplate was monitored for growth using an automated microplate reader at 24-h intervals for 3 days. Isolates that demonstrated growth at the following minimum metal ion concentrations were classified as resistant using previously described criteria (Dziewit et al., 2015): (i) 20 mM V⁵⁺, (ii) 1 mM As³⁺, Cd²⁺, Co²⁺, Cu²⁺, Ni²⁺, Zn²⁺ or Cr⁶⁺, and (iii) 0.1 mM Hg²⁺.

Plasmid DNA was isolated using the alkaline lysis procedure (Birnboim and Doly, 1979) or purified by CsCl-ethidium bromide gradient centrifugation (Sambrook and Russell, 2001). The visualization of mega-sized replicons was achieved by in-gel lysis and DNA electrophoresis (Eckhardt, 1978; Hynes and McGregor, 1990). Plasmid DNA was also isolated using a Plasmid Mini kit (A&A Biotechnology), GeneMATRIX Miniprep DNA Purification Kit (EUR_x) and GeneJET Plasmid Miniprep Kit (Thermo Fisher Scientific). Routine DNA manipulation was performed using standard methods (Sambrook and Russell, 2001). DNA amplification by PCR was performed in a Mastercycler (Eppendorf) using synthetic oligonucleotides (Supplementary Table S1), High Fidelity Taq DNA polymerase (Qiagen) or Phusion Hi-Fidelity DNA polymerase (Thermo Fisher Scientific), dNTPs and appropriate template DNAs, as described previously (Bartosik et al., 2003).

Chemical transformation of E. coli cells was performed by a standard method (Kushner, 1978). Plasmid DNA was introduced into Paracoccus spp. strains by triparental mating using helper E. coli strain DH5α carrying plasmid pRK2013 (containing the transfer system of plasmid RP4) (Ditta et al., 1980), as described previously (Lasek et al., 2018).

REP systems of extrachromosomal replicons of the analyzed strains were cloned in vector pABW1 (Bartosik et al., 1997) which cannot replicate in Paracoccus strains or other tested recipients (Supplementary Table S1). These shuttle plasmid constructs were introduced into the following strains: (i) class Alphaproteobacteria—Agrobacterium tumefaciens LBA 288R (Hooykaas et al., 1980), Paracoccus aminophilus JCM 7686R (without plasmids pAMI2, pAMI3) (Dziewit et al., 2014), Paracoccus aminovorans JCM 7686R (Czarnecki et al., 2017), Paracoccus pantotrophus KL100 (Jordan et al., 1997), (ii) class Betaproteobacteria—Achromobacter sp. LM16R (Dziewit et al., 2015) and (iii) class Gammaproteobacteria—Cronobacter sakazakii ATCC 29544 (Iversen et al., 2008) and Pseudomonas sp. LM6R (Dziewit et al., 2015).

Trap plasmids pMAT1 and pMEC1 (Bartosik et al., 2003; Szuplewska and Bartosik, 2009; Gay et al., 1985; Schneider et al., 2000) were used for the identification of functional TEs of P. yeei (Supplementary Table S1). These plasmids were transferred from E. coli TG1 to P. yeei strains by triparental mating. The “capture” of TEs was verified by PCR using primers specific to the selection cassettes, with DNA isolated from Tc^r and Suc^r mutants as the template (Supplementary Table S1). A total of 1,304 clones were tested this way. The amplified DNA fragments were subjected to electrophoretic analysis, which permitted the detection of DNA elements (TEs), ranging in size from ~0.8 to 1.5 kb, embedded in the cassettes. DNA sequencing of the termini of individual TEs was used to identify the complete elements within the genome sequences of their parental strains. All trapped TEs were sequenced and compared with the ISfinder database (Siguier et al., 2006).

To confirm that the identified URE modules are responsible for ureolytic activity, mutational analysis was performed using the gene replacement method. Mutations (insertion of a Km^r cassette) were introduced into URE type 1 modules, either in the 5′ region of the Ureα subunit gene or within the nikR gene. Both mutations were created in the strain CCUG 32053 (which contains only URE type 1) and the strain CCUG 13493 (which also contains a URE type 2 module). The following procedure was used to obtain ureα::Km^r and nikR::Km^r mutants: (i) DNA fragments containing ureα-ureβ-ureγ or part of the region encoding NikR and nickel transporters were amplified by PCR using CCUG 32053 DNA as the template, (ii) these DNA fragments were cloned into the SacI and SphI sites of the mobilizable suicide vector pDS132 (Philippe et al., 2004), (iii) the Km^r cassette from plasmid pDIY-Km (Dziewit et al., 2011) was inserted into the AleI site within these fragments to obtain constructs pDS132-ureα and pDS132-nikR (Supplementary Table S1); (iv) these plasmids were then introduced into strains CCUG 32053R and CCUG 13493R by conjugation. Mutant clones ureα::Km^r and nikR::Km^r were selected on medium supplemented with kanamycin (for recombinant selection) and sucrose (for counter-selection). The introduced mutations were confirmed by DNA sequencing. Ureolytic activity of the obtained mutants was tested following growth on minimal medium supplemented with 1% urea, using rapid urease tests. Complementation of the mutants was performed using plasmids pBBR-ureα or pBBR-nikR, which were prepared by cloning PCR fragments containing ureα-ureβ-ureγ or part of the region encoding NikR and nickel transporters into vector pBBR-MCS5 (Gm^r) (Kovach et al., 1995) (Supplementary Table S1).

The P. yeei genomes were sequenced using a combination of Oxford Nanopore and Illumina technologies, as described previously (Lasek et al., 2018).

Sequence annotation and bioinformatic analyses (identification of tRNA genes and rRNA operons, relaxases—MOB, toxin-antitoxin modules, transposable elements—TE) were performed as described previously (Lasek et al., 2018). Potential ICE and IME elements were predicted using ICEfinder (Liu et al., 2019).

Categories for clusters of orthologous groups (COGs) were assigned to each protein by performing a local RPS-BLAST search against the COG database (last updated on January 22, 2015). A threshold e-value of 1e-5 was used, and only the top BLAST hits were taken into account (Tatusov et al., 2003).

The core genome of Paracoccus spp. and P. yeei species was defined based on the complete genomic sequences of P. aminophilus JCM 7686 (Dziewit et al., 2014), P. aminovorans JCM 7685 (Czarnecki et al., 2017), P. contaminans RKI16-01929T (Aurass et al., 2017) and P. denitrificans PD1222 (NC_008686-8), and eleven P. yeei strains: CCUG 32053 (Lasek et al., 2018), FDAARGOS_252 (NZ_CP020440–47), FDAARGOS_643 (NZ_CP044078-82), TT13 (Lim et al., 2018) and CCUG 13493, CCUG 17731, CCUG 32052, CCUG 32054, CCUG 46822, CCUG 54214, LM20 (Dziewit et al., 2015). Proteins encoded within the genomes were used in all-against-all BLASTp searches utilized by a stand-alone version of OrthoVenn2 (Xu et al., 2019), using an e-value of 1e-15 as the threshold and inflation level of 1.5. Based on this analysis, all proteins were clustered into groups reflecting similarity and designated as core proteins when encoded by all eleven genomes or as singletons when the respective gene was identified in only a single genome. A similar approach with a different manner of grouping was applied to identify strain-specific gene clusters.

Protein datasets from the Virulence Factor Database VFDB (Chen et al., 2016), the Comprehensive Antibiotic Resistance Database CARD (Alcock et al., 2020) and Virulence Factors Database VICTORS (Sayers et al., 2019), all downloaded on April 1 2021, were used for BLASTp searches (e-value threshold of 1e-30 and 80% of query coverage per HSP) with the proteomes of P. yeei and four other Paracoccus spp. strains listed above to identify putative P. yeei species virulence determinants.

EasyFig (Sullivan et al., 2011) was used to perform comparative genomic analyses (including genomes and plasmidomes comparisons) and visualize the results. RNA secondary structures were predicted by in silico folding using Mfold software (Zuker, 2003).

Phylogenetic analysis of the genus Paracoccus was based on alignment of concatenated nucleotide sequences of selected core genes: atpD, dnaA, dnaK, gyrB, recA, rpoB, and thrC (homologous genes of Roseobacter denitrificans OCh 114 were used as an outgroup), as described previously (Lasek et al., 2018). Nucleotide alignments were obtained with Mafft (Katoh et al., 2019). Then, concatenated genes were analyzed with ModelTest-NG (Darriba et al., 2020), checking all models to select the best-fit nucleotide substitution model. The selected substitution model was applied in RaxML-NG (Kozlov et al., 2019) with 2,000 regular bootstrap replicates performed on the best Maximum Likelihood (ML) tree selected from 100 independently generated ML starting trees.

The nucleotide sequences of the P. yeei chromosomes and extrachromosomal replicons were deposited in GenBank (NCBI), with the following accession numbers: (i) CCUG 13493—CP038080 (chromosome), CP038073 (plasmid pYEE13493P1), CP038074 (pYEE13493P2), CP038075 (pYEE13493P3), CP038076 (pYEE13493P4), CP038077 (pYEE13493P5), CP038078 (pYEE13493P7), CP038079 (pYEE13493P8), CP038081 (pYEE13493P6); (ii) CCUG 17731—CP038042 (chromosome), CP038035 (pYEE17731P1), CP038036 (pYEE17731P2), CP038037 (pYEE17731P3), CP038038 (pYEE17731P4), CP038039 (pYEE17731P5), CP038040 (pYEE17731P6), CP038041 (pYEE17731P7); (iii) CCUG 32052—CP038090 (chromosome), CP038082 (pYEE32052P1), CP038083 (pYEE32052P2), CP038084 (pYEE32052P3), CP038085 (pYEE32052P4), CP038086 (pYEE32052P5), CP038087 (pYEE32052P6), CP038088 (pYEE32052P7), CP038089 (pYEE32052P8); (iv) CCUG 32054—CP038095 (chromosome), CP038091 (pYEE32054P6), CP038092 (pYEE32054P7), CP038093 (pYEE32054P8), CP038094 (pYEE32054P9), CP038096 (pYEE32054P1), CP038097 (pYEE32054P2), CP038098 (pYEE32054P3), CP038099 (pYEE32054P4), CP038100 (pYEE32054P5); (v) CCUG 46822—CP038056 (chromosome), CP038043 (pYEE46822P1), CP038044 (pYEE46822P2), CP038045 (pYEE46822P3), CP038046 (pYEE46822P4), CP038047 (pYEE46822P5), CP038048 (pYEE46822P6), CP038049 (pYEE46822P7), CP038050 (pYEE46822P8), CP038051 (pYEE46822P9), CP038052 (pYEE46822P10), CP038053 (pYEE46822P11), CP038054 (pYEE46822P12), CP038055 (pYEE46822P13); (vi) CCUG 54214—CP038061 (chromosome), CP038057 (pYEE51214P1), CP038058 (pYEE51214P2), CP038059 (pYEE51214P3), CP038060 (pYEE51214P4); (vii) LM20—CP038072 (chromosome), CP038062 (pLM20P1), CP038063 (pLM20P2), CP038064 (pLM20P3), CP038065 (pLM20P4), CP038066 (pLM20P5), CP038067 (pLM20P6), CP038068 (pLM20P7), CP038069 (pLM20P8), CP038070 (pLM20P9), CP038071 (pLM20P10). The nucleotide sequences of newly identified ISs (ISPye73—ISPye80) and MITEPye1 were deposited in the ISfinder database (Siguier et al., 2006).

Seven (five clinical and two environmental) isolates of P. yeei were selected for detailed characterization (Table 1). Clinical strains were sourced from two continents—they were isolated in different geographical locations, at different times and from different clinical cases (Table 1). The other strains were isolated either from an anthropogenic environment (CCUG 54214; a metal surface), or from the natural environment (LM20; organic-rich black shale from the Lubin copper mine in Poland; Dziewit et al., 2015). The clinical isolate CCUG 32053, characterized in a previous study (Lasek et al., 2018), was used as a reference strain.

The genomes of the seven strains of P. yeei were fully sequenced and characterized. Their genome sizes range from 4,423,927 to 4,829,807 bp (average 4,610,003 bp). All are multipartite, containing from 4 (CCUG 54214) to 13 (CCUG 46822) extrachromosomal replicons (ECRs), that range in size from 3.5 kb (pYEE46822P1 of CCUG 46822) to 485 kb (pYEE17731P7 of CCUG 17731) (Table 3; Figure 1B; Supplementary Figure S1). In total, the strains contain 59 ECRs (11,665,529 bp), which on average constitute ~25% of the size of each genome (Table 3; Figure 1B). The physicochemical parameters of the genome sequences were determined, including the %GC content (average 67%), the number of protein-coding sequences (CDSs) (average 4,200) and pseudogenes (average 190). The GC content of the ECRs ranges from 53% (pYEE46822P4, CCUG 46822) to 69% (pYEE13493P5, CCUG 13493). The number of CDSs present in ECRs ranges from 849 (CCUG 54214), accounting for 20% of the entire genome, to 1,270 (CCUG 17731) (28%) (Supplementary Table S3). The rRNA operons and most of the tRNA genes are localized within the chromosomes, although some additional tRNA genes (2–3) were identified within ECRs (Table 3).

The analyzed P. yeei strains contain numerous ECRs (59 in total)—from 4 (CCUG 54214) to 13 (CCUG 46822) (Figure 1B; Supplementary Figure S1). These are mainly large DNA molecules, constituting a significant proportion of their host strain genomes—from 22% (CCUG 54214) to >28% (CCUG 13493) (Figure 1B). The replication system (REP) and the gene encoding the replication initiation protein (Rep) within each ECR were identified. The REP modules (and entire plasmids) were classified based on the presence of specific amino acid (aa) sequence signatures in the Rep proteins (Figure 2A). Replicons of the DnaA-like, RepB (RepB1 and RepB2), RepABC, RepC, Rep_3 and HTH_36 families were identified, all of which are common among Paracoccus spp. and other Alphaproteobacteria (Figures 1C, 2).

The ECRs carrying REPs of a given type are generally similar in size, with the exception of RepB1 and RepB2 replicons, in which greater size variation is observed (Figure 2B). The REPs of the DnaA-like, RepB and RepABC types (occurring in each P. yeei strain) are characteristic for large replicons, which is consistent with previous observations (e.g., Lasek et al., 2018). Detailed characterization of individual ECRs, including their size, nucleotide sequence GC content, coding capacity and tRNA gene distribution, is presented in Supplementary Table S3.

To test the host range of individual ECRs representing the different replicon families, their REP regions were cloned in mobilizable shuttle plasmids and introduced by triparental mating into recipient strains from different classes of Proteobacteria (listed in Section 2). This analysis revealed that all of the analyzed REPs function solely in Alphaproteobacteria, so have a relatively narrow host range.

The P. yeei ECRs contain numerous loci enabling their stable maintenance in bacterial cells and populations. The DnaA-like, RepB, and RepABC replicons all carry type I partitioning systems, encoding Walker-type ATPases (Thomas, 1981) (Figure 2A). The vast majority of ECRs also have class II toxin-antitoxin (TA) systems, mainly encoding toxins of the VapC family (18 toxins), as well as toxins of the RelE/ParE (14), PhD/YefM (14) or parDE (8) families. None of the ECRs encode a complete type IV secretion system (T4SS), suggesting that these replicons are not self-transmissible. However, 16 of them contain genes encoding predicted relaxases of the Mob_C (14), Mob_Q (1), and Mob_HEN (1) families (Garcillán-Barcia et al., 2009), which are typical components of genetic modules, enabling plasmid mobilization for conjugative transfer (MOB).

Further analysis was conducted to verify whether the P. yeei genomes contain ECRs that meet the criteria for chromids (secondary chromosomes)—essential replicons of plasmid origin. This revealed that all the DnaA-like, RepB1 and RepB2 replicons can be considered chromids (Figure 1C), since they (i) carry plasmid-type REP modules, (ii) have a nucleotide sequence composition that is close to that of the chromosome (< 1% difference in GC content), and (iii) possess a set of core genes (in one copy in the genome) characteristic for the entire genus (Harrison et al., 2010). The genes carried by these chromids participate in several important metabolic pathways, e.g., encoding acetyl-CoA C-acyltransferase, and proteins likely to be involved in biosynthesis of the lipoyl cofactor, in oxidative metabolism, leucine metabolism, and electron transport via the respiratory chain (complex I).

The Rep proteins of the DnaA-like, RepB1 and RepB2 chromids display a high level of aa sequence identity within a given replicon group. These replicons are also more highly conserved in terms of their structure and genetic load (Supplementary Table S5). In contrast, the replication initiators of the RepABC, Rep_3 and HTH_36 families (and their entire replicons) show much greater sequence variability. Many of these replicons appear to be unique among P. yeei strains, suggesting their recent horizontal acquisition (Supplementary Table S5).

Phylogenetic analysis showed that the two environmental isolates (CCUG 54214 and LM20) are situated on different branches of the P. yeei phylogenetic tree and are closely related to the clinical isolates. Strains LM20 and CCUG 46822 form a separate cluster, while the other five strains can be divided into two phylogenetic groups: (i) CCUG 54214 and CCUG 13493, and (ii) CCUG 17731, CCUG 32052 and CCUG 32054 (Figure 1A).

Comparative genomic analysis of the P. yeei strains (including clinical isolates CCUG 32053, FDAARGOS_643, FDAARGOS_252, and TT13) revealed a conserved genome structure, preserved ECR composition and synteny of gene arrangement in this species (Supplementary Table S5). There is a particularly high degree of synteny between the genomes of the strains LM20 and CCUG 46822 (Figure 3). The ECRs of these strains are nearly identical, although they display minor rearrangements (Supplementary Table S5).

Our analysis identified two chromosome structure variants (one found in strains CCUG 17731, CCUG 32052, CCUG 32053, CCUG 32054, CCUG 54214, FDAARGOS_643, TT13, and the other in CCUG 13493, CCUG 46822, LM20, FDAARGOS_252), which correlates with the phylogeny of the strains (Figure 3; Supplementary Table S5). The differences result mainly from inversions or translocations of fairly large DNA segments (Figure 3).

The DnaA-like replicons show a high degree of conservation, although with some rearrangements, e.g., translocation of DNA segments in strains CCUG 32053 and CCUG 13493, sometimes with inversion of the transferred DNA segment (Supplementary Table S5).

Analysis of the proteome of P. yeei determined in silico permitted the assignment of individual proteins to functional COG groups. For all genomes a similar proportion of the encoded proteins belong to each category, indicating a conserved genome structure and limited variability within the species. Interestingly, significant similarity in the distribution of COG groups was found between the chromosomes and ECRs, suggesting that ECRs may contain a large amount of genetic information of chromosomal origin (Figure 4).

In all cases, a considerable number of genes are of unknown function (R and S categories) or could not be classified into any COG group: 27–31% for chromosomes and 27–32% for ECRs. The largest group of classified proteins comprises those related to amino acid transport and metabolism (category E; 9–10%), followed by carbohydrate metabolism (G; 5–6%), transcription (K; 6%), transport and metabolism of inorganic ions (P; 5–6%), energy acquisition processes (C; 5%), and the synthesis of cell wall, membrane, and envelope components, including capsules (M; 5%) (Supplementary Figure S2). In the case of ECRs, the largest fraction of classified proteins was assigned to categories E (7–10%), G (8–11%), P (6–8%), and K (6–8%) (Figure 4B). There are no significant differences between clinical and environmental isolates in the distribution of COG categories.

Transport-related genes constitute a dominant group of functional genes in the P. yeei genomes—from 530 (CCUG 32052) to over 600 CDSs (CCUG 17731). These genes encode various types of transporters: ABC, MSF (Major Facilitator Superfamily), P-type, RND (Resistance-Nodulation-Division), TRAP (Tripartite ATP-independent Periplasmic), and TTT (Tripartite Tricarboxylate Transporters) (Table 4; Supplementary Table S6). ECRs (mainly DnaA-like and RepB chromids) contain over 35% of the ABC transporter genes. Less common in ECRs are MFS, RND and TTT transporter genes. Only a few genes encoding proteins with similarity to P-type ATPase integral membrane transporters were found in the analyzed genomes.

The ECRs also carry numerous genes that may increase the adaptability of their host strains in the environment. All types of chromids and RepABC plasmids encode potential iron or carbohydrate transporters. Siderophore receptor genes were identified in all RepB2 plasmids, including the FhuADCB ferric hydroxamate transporter, which is essential for the uptake of Fe³⁺-aerobactin. All RepB1 replicons carry a mauMGJCBDEAF gene cluster involved in the utilization of methylamine (Baker et al., 1998). Organisms that utilize C1 compounds as their sole carbon and energy sources are described as methylotrophs and this trait is common in bacteria of the genus Paracoccus (Czarnecki and Bartosik, 2019). Methylotrophs play an important role in global carbon, nitrogen and sulfur cycling, and they have been successfully employed in the bioremediation of contaminated soils (Chistoserdova, 2011).

The RepB1 replicons of three phylogenetically closely related P. yeei strains (CCUG 17731, CCUG 32052, CCUG 32054) (Figure 1A) contain a putative operon (phnCDEFGHIJKLMNOP) encoding carbon-phosphorus lyase, which permits the utilization of phosphorus from a wide range of stable phosphonate compounds containing a C-P bond (Makino et al., 1991; Horsman and Zechel, 2017; Amstrup et al., 2023). In numerous environments, inorganic phosphate, a vital nutrient, is available in very limited amounts, compelling microorganisms to rely on alternative phosphorus sources for survival (Podzelinska et al., 2009). Therefore, the ability to uptake and breakdown phosphonates would be advantageous.

Enzymes catalyzing rhamnose biosynthesis (RfbABCD) are encoded by all RepB2 replicons (CCUG 32052 has an additional copy of the rfbABCD operon in Rep_3 plasmid pYEE32052P3). Rhamnose is present in the O-antigens of many gram-negative bacteria, forming part of the lipopolysaccharide (LPS). This deoxy-hexose sugar is also present in capsular polysaccharides, covalently bound to the cell wall, and in exopolysaccharides that are loosely associated with the cell wall (Marolda and Valvano, 1995). The P. yeei strains analyzed in this study display characteristic mucoid growth, probably related to the presence of a polysaccharide capsule (data not shown).

The ECRs also carry genes for proteins involved in (i) the cobalamin (vitamin B12) biosynthetic pathway (cobAWNGHIJKLMBF; DnaA-like chromids), (ii) ectoine biosynthesis (ehuABCD; RepB1 chromids), (iii) propionate catabolism (mmgE/prpD; several RepB1, RepB2, and RepABC replicons), (iv) disulfide bond formation (dsbA, dsbB, and dsbE), and (v) resistance to arsenic compounds (permease Arc3, arsenate reductase ArsC, flavoprotein ArsH, repressor ArsR, and transporter ArsJ). Notably, all RepB1 chromids carry a complete conserved set of crt genes (crtXYIBZE-idi) for carotenoid synthesis (Maj et al., 2013). However, colonies of the strains with these chromids lack the characteristic color associated with carotenoid pigment production and these compounds could not be detected by HPLC analysis (data not shown).

In summary, the DnaA-like, RepB1 and RepB2 chromids and the RepABC plasmids carry the largest number of genes of adaptive potential. The genetic load of smaller replicons (HTH_36 or Rep_3) consists mainly of genes of unknown function. The distribution of selected genes in individual ECRs is presented in Supplementary Table S3.

A search for transposable elements (TEs) in the P. yeei genomes revealed the presence of over 900 tnp genes (complete or partial) encoding transposases from 18 families of insertion sequences (ISs): IS3, IS5, IS6, IS21, IS30, IS66, IS91, IS110, IS256, IS481, IS630, IS701, IS1182, IS1202, IS1380, IS1595, ISL3, and ISNCY (Figure 5; Supplementary Table S7). The number of tnp genes in individual genomes ranges from 84 (CCUG 32054) to 150 (CCUG 32052 and LM20). The dominant elements are representatives of the IS3, IS5, IS110, and IS256 families. Several IS families were identified in only one or a few isolates, e.g., IS6—CCUG 17731 and CCUG 32052, IS91—CCUG 46822, IS630—CCUG 32054 and LM20, IS1202—CCUG 46822, IS1380—CCUG 46822, IS1595—CCUG 13493, CCUG 46822 and LM20, and ISNCY—CCUG 13493 and CCUG 32054 (Supplementary Table S7).

Numerous ISs are localized within ECRs. Surprisingly, most of the largest ECRs (DnaA-like and RepB chromids), except for those from CCUG 46822 and LM20, contain relatively few ISs, indicating a more conserved replicon structure and the presence of essential genes. The replicons richest in these elements are RepABC plasmids, which in most cases carry almost 50% of the tnp genes of a given strain (Supplementary Table S7).

In the chromosomes and ECRs of strains CCUG 46822 and LM20 several repetitive IS arrays were identified with a structure resembling composite transposon (Supplementary Table S7). These DNA regions are flanked by divergently oriented complete isoforms of ISPye46 (1361 bp; IS110 family) and contained within them are tnp genes of the IS256 family (related to ISPye43) and/or the IS5 family (related to ISPye12). Each copy of ISPye46 is flanked by identical DRs (5′-TG-3′), so it is not clear whether these IS110-IS5/IS256-IS110 arrays constitute functional TEs or are the result of preferential transposition of ISPye46 into ISs representing the IS256 or IS5 families. Interestingly, ISPye46 is dominant in the analyzed P. yeei genomes (47 copies in total); however, most of the copies are present in strains CCUG 46822 (19 copies) and LM20 (24 copies) (Figure 5; Supplementary Table S7). No related IS110-IS5/IS256-IS110 arrays were identified in the genomes of other bacteria by BLASTn analysis.

Interestingly, a 159-bp-long sequence element was identified in most chromosomes of the analyzed P. yeei strains (Supplementary Table S7). These elements show high sequence similarity (96–100% identity), lack any open reading frames, occur in intergenic regions at different locations, and are flanked by short (4 bp) direct repeats (DRs), suggesting their acquisition through transposition (Supplementary Table S7). Additionally, they contain inverted repeats at both termini (TIRs; IRL and IRR) (18 bp) with significant similarity to the TIRs of insertion sequences of the IS1182 family, including ISPye18, previously identified in the CCUG 32053 genome (Lasek et al., 2018) (Figure 6A). Analysis of the distribution of ISPye18 in P. yeei genomes revealed the presence of intact isoforms of this element in strains FDAARGOS_252 (1 copy; plasmid 4) and CCUG 32053 (1 copy; plasmid pYEE3). However, related ISs from the IS1182 family were also identified in the genomes of CCUG 13493, CCUG 46822, CCUG 54214, and LM20 (Supplementary Figure S3).

Considering (i) the potential mobility of these DNA regions, (ii) their small size, (iii) the presence of TIRs, (iv) DRs, and (v) the absence of a transposase gene, it seems likely that they constitute a novel group of MITE-type non-autonomous TEs (Delihas, 2011), that are most probably mobilized for transposition by ISs of the IS1182 family. This element was named MITEPye1 (ISfinder).

The activity of TEs of P. yeei was investigated in vivo using trap plasmids that allow positive selection of transposition events. Two trap plasmids were used, pMEC1 and pMAT1, containing the cI-tetA and sacB selection cassettes, respectively. Transposition of TEs into these cassettes resulted in the appearance of new phenotypes—the cells carrying a mutated plasmid become resistant to tetracycline (pMEC1) or to sucrose (pMAT1) (see Section 2 for details).

Despite the presence of a very large number of transposase genes in the P. yeei genomes, use of the entrapment vectors led to the identification of only six distinct elements: (i) three ISs previously identified in silico in the CCUG 32053 genome—ISPye2 (IS5 family, IS5 group), ISPye38 (IS5 family) and ISPye41 (IS5 family, IS903 group), (ii) an isoform of IS1247 of Xanthobacter autotrophicus (IS1380 family), and (iii) two novel elements, designated ISPye79 and ISPye80, representing the IS5 and IS427 groups within the IS5 family, respectively.

The distribution of every defined IS of P. yeei (ISfinder and this study) across all P. yeei genomes is shown in Supplementary Figure S3. The most abundant elements are (i) ISPye46−47 copies in 5 strains (most copies present in CCUG 46822 and LM20), (ii) ISPye53−45 copies in seven strains, and (iii) ISPye41−30 copies in six strains (most copies present in CCUG 17731 and CCUG 32053). Many of the ISs are present in only 1 copy, e.g., (i) ISPye13 (chromosome FDAARGOS_252), (ii) ISPye63 (chromosome of CCUG 32052), and (iii) IS1247a (chromosome of CCUG 46822).

The core genomes of Paracoccus spp. and of P. yeei were defined based on the complete genomic sequences of P. aminophilus JCM 7686 (Dziewit et al., 2014), P. aminovorans JCM 7685 (Czarnecki et al., 2017), P. contaminans RKI16-01929T (Aurass et al., 2017), P. denitrificans PD1222, plus eleven P. yeei strains. Comparative genomic analysis was used to define (i) the core genome of bacteria of the genus Paracoccus spp., (ii) genes specific to P. yeei (not found in other Paracoccus spp. with complete genomes), and (iii) genes that are unique to individual P. yeei strains (singletons). The collected data, also showing the distribution of these genes in the individual replicons of each P. yeei strain, are presented in Supplementary Table S8.

This analysis revealed that the core genome of Paracoccus spp. consists of 1,647 protein-encoding genes, of which 116 are located within ECRs of P. yeei. The total number of P. yeei-specific genes is 250 and the number of singletons ranges from 138 (CCUG 46822) to 332 (FDAARGOS 252) (Supplementary Table S8). In some strains, the majority of singletons occur in the chromosome, while in others (CCUG 13493, CCUG 17731) they are predominantly in ECRs. The vast majority of singletons encode proteins of unknown function, although a predicted role in adaptation could be assigned to some genes from the P. yeei environmental isolates. The chromosome of LM20 encodes (i) the membrane protein YeeE, which mediates the uptake of thiosulfate as an inorganic sulfur source for cysteine synthesis (Tanaka et al., 2020), (ii) mercury transporter Mer, and (iii) the sensor kinase KdpD, which senses potassium levels, and the kdpABC operon, which encodes a high-affinity potassium uptake system crucial for bacterial survival in low potassium environments (Laermann et al., 2013; Ali et al., 2017). In the CCUG 54214 genome, many transport-related genes (ABC, RND, TTT types) were distinguished among the singletons, as well as genes for toxin-antitoxin systems (e.g., HipA-, PrlF-, Yha- family), similarly to clinical strains.

The largest group of P. yeei-specific genes (found in chromosomes and DnaA-like, RepB1, RepB2 and RepABC replicons) encode transposases from different IS families. The second largest group of genes encodes transporters, mostly of the ABC type, which are located in chromosomes and in DnaA-like and RepB1 chromids. Other abundant genes, found mainly in RepB1 chromids, encode transcription regulators of the LuxR, CopG, TetR/AcrR, MarR, LacI, AraC, TetR, DeoR/GlpR, and NikR families (Supplementary Table S8).

Examples of proteins encoded by other P. yeei-specific genes include a urease enzyme, components of the trehalose, rhamnose, and cobalamin biosynthetic pathways, chaperone proteins, and factors possibly involved in the cell's response to stress conditions (Supplementary Table S8).

The proteomes of a number of Paracoccus strains were analyzed in order to identify potential virulence factors: (i) the seven P. yeei strains described in this study, (ii) other P. yeei strains whose complete genomic sequences are available (CCUG 32053, TT13, FDAARGOS_262, FDAARGOS_643), as well as (iii) several environmental isolates of the genus Paracoccus, including the type strain P. denitrificans PD1222. These protein sets were compared against databases containing bacterial virulence factors and antibiotic resistance genes (see Section 2 for details). Relatively few of the identified determinants were unique to P. yeei (all clinical isolates). The pool of identified genes includes those encoding six types of protein, that were previously mentioned by Lasek et al. (2018), i.e., (i) methionine sulfoxide reductase (msrA1, msrA2, msrB), (ii) diguanylate cyclase, (iii) superoxide dismutase, (iv) sugar transferase, (v) type IV secretion system components, and (vi) urease (Supplementary Tables S3, S9).

The chromosomally-located msr genes of CCUG 32052 and CCUG 46822 encode predicted proteins sharing 61–100% aa sequence similarity, while the homologous proteins encoded by other P. yeei strains are less well-conserved (38–55% aa identity). Msr proteins participate in the detoxification of reactive oxygen intermediates (both prokaryotic and eukaryotic cells lacking Msr are sensitive to oxidative stress) (Denkel et al., 2011). Moreover, MsrA has been shown to be critical for the survival of Erwinia chrysanthemi (El Hassouni et al., 1999), Mycoplasma genitalium (Das et al., 2012), and Helicobacter pylori (Alamuri and Maier, 2006) in their infected hosts. Genes encoding components of a putative type IV secretion system (T4SS) were identified only in CCUG 17731 (30–50% aa identity) and in CCUG 13493 (25–50% aa identity)—within RepABC plasmids in both cases. However, these were not complete modules, unlike the one previously identified in CCUG 32053 (Lasek et al., 2018).

The remaining virulence factors, present in all analyzed genomes, include the gene encoding diguanylate cyclase, an enzyme synthesizing cyclic diguanylate (c-di-GMP)—an important bacterial second messenger. This signaling molecule is involved in the regulation of a number of complex physiological processes, including biofilm formation and motility, which affect the pathogenesis of many bacteria (Tamayo et al., 2007).

Genes encoding several other putative virulence factors were identified within ECRs, although they were not unique to P. yeei. These include proteins participating in the following: (i) numerous ABC and RND transport systems (DnaA-like, RepB, RepABC replicons), (ii) siderophore synthesis (RepB2 replicons, chromosomes), (iii) polysaccharide capsule production (RepB2 replicons), and (iv) trehalose synthesis (DnaA-like type replicons).

All P. yeei strains contain a set of genes (URE) involved in the synthesis of the metalloenzyme urease, which hydrolyzes urea to ammonia and CO₂. These genes are located within essential replicons—DnaA-like chromids. They display synteny and a high degree of nucleotide sequence identity (98%), indicating their biological importance. The role of urease in bacterial pathogenesis is well-documented (e.g., Mobley, 1996; Dupuy et al., 1997; Burne and Chen, 2000; Berutti et al., 2014; Rutherford, 2014; Graham and Miftahussurur, 2018; Zhou et al., 2019; Minami et al., 2021).

Besides genes for the main urease subunits (UreA—γ subunit, UreB—β subunit, and UreC—α subunit) and accessory proteins (ureD, ureE, ureF, and ureG), the URE gene cluster also contains the nikR gene (encoding a predicted nickel-dependent transcription factor) and numerous ABC transporter genes, likely to be involved in the uptake of Ni²⁺ (an essential cofactor for urease; Rutherford, 2014). It is probable that NikR regulates the expression of the transporter genes (Figure 7). Interestingly, in P. yeei strain CCUG 32052, part of the URE gene cluster, containing the ureABC and ureD genes, has been duplicated and transferred to another ECR (RepABC; pYEE32052P4) (99% nucleotide sequence identity to the corresponding URE genes) (Figure 7).

The URE gene cluster of DnaA-like replicons (designated URE type 1) (Figure 7) is unique to and characteristic of P. yeei strains. However, in two strains (CCUG 13493 and CCUG 32052) another set of ure genes (URE type 2; Figure 7), possibly associated with urease synthesis, was identified. These genes are present in RepABC plasmids in both strains (pYEE13493P4 and pYEE32052P5, respectively) and their genetic organization is different from that of URE type 1. In the URE type 2 cluster, the γ and β urease subunits are encoded by a single fusion gene ureBC, the order and transcriptional orientation of the accessory genes ureFEDG is different, there are fewer transporter genes and the nikR gene is absent (Figure 7).

URE type 2 represents a novel type of ureolytic gene cluster which has a different evolutionary origin from URE type 1. The two clusters share some sequence similarity between the genes encoding the main urease subunits, but no similarity was found in case of the accessory genes and their arrangement is different.

The ureolytic activity of three P. yeei strains (CCUG 32053, CCUG 32052—containing the duplicated part of URE type 1, and CCUG 32052—also carrying URE type 2; Figure 7) was examined using Christensen's differential medium without (control) and with the addition of urea (1%) (see Section 2 for details). All the strains displayed ureolytic activity, which was induced by the presence of urea in the medium.

To unequivocally confirm the role of the type 1 URE gene cluster in the observed ureolytic activity, mutations were introduced by insertion of a kanamycin resistance cassette into the ureC (encoding the α subunit of urease) or nikR genes of CCUG 32053 and CCUG 13493. These mutations completely abolished ureolytic activity in these strains, indicating the importance of both ureC and nikR in determining this phenotype. The lack of ureolytic activity in mutated strain CCUG 13493, which also contains URE type 2, indicates that this second cluster does not produce an active urease under the tested conditions and its genes are unable to complement mutations introduced in URE type 1.