Five copper homeostasis gene clusters encode the Cu-efflux resistome of the highly copper-tolerant Methylorubrum extorquens AM1

Relevant characteristics of PAFD and PPFM

The strains analyzed in this work (Table 1) are model PAFD and PPFM, presumed to come from niches not polluted with metals, widely used in laboratory and field studies. All the strains were kindly donated by Professor Esperanza Martínez Romero (Centro de Ciencias Genómicas, UNAM, México City, México), and were carefully maintained in 20% glycerol stored at −70 °C. All the PAFD were isolated from the inside of root nodules of different species of leguminous plants growing in a variety of habitats in different countries (Dall’Agnol et al., 2013; Galibert et al., 2001; Herrera-Cervera et al., 2006; Johnston & Beringer, 1974; Martínez-Romero et al., 1991; Saeki & Kouchi, 2000; Segovia, Young & Martínez-Romero, 1993; Trlnlck, 1980). The studied PPFM are cosmopolitan freestyle bacteria isolated from different environments around the world. This group includes two Mr. extorquens strains and three Mb. species. Mb. nodulans ORS2060^T is a facultative methylotroph and diazotroph isolated from the inside of root nodules of Crotalaria podocarpa in Senegal. Although it is a non-pigmented species, to simplify the results, it was included in the PPFM group (Ito & Iizuka, 1971; Jourand et al., 2004; Patt, Cole & Hanson, 1976; Peel & Quayle, 1961; Urakamit & Komagata, 1984). The ability of these bacteria to grow from methanol as the sole carbon source and the pink-pigmented phenotype were verified in this study. The sequenced genomes of PAFD and PPFM can be freely accessed through NCBI genome resources (https://www.ncbi.nlm.nih.gov/home/genomes/). A list of accession numbers is included in Table S1.

Media cultures

Bacterial cells maintained at −70 °C were propagated in peptone-yeast (PY)-agar medium containing 5 g/L peptone, 3 g/L yeast extract, and 15 g/L agar. After sterilization, 10 ml/L 0.7 M CaCl₂ was added.

Minimal medium (Mm) is a chemically defined medium prepared from three solutions (A, B, and C) and sterilized separately. Solution A contained 1.620 g/L sodium succinate hexahydrate as a carbon source, 0.534 g/L NH₄Cl as a nitrogen source, 0.219 g/L K₂HPO₄, and 0.1 g/L MgSO₄, and its pH was adjusted to 6.8 before sterilization. Agar (15 g/L) was added and then sterilized in an autoclave. Solution B contained filter-sterilized 0.025 g/5 ml FeCl₃·6H₂O, and solution C contained 0.7 M CaCl₂·2H₂O (autoclaved). One milliliter of B solution and two milliliters of C solution were added to one liter of previously sterilized A solution.

Estimation of minimal inhibitory concentrations (MICs) of CuCl₂

We followed the definition and standardized dilution method set by the British Society of Antimicrobial Chemotherapy and updated by Andrews (2001), and by the European Committee for Antimicrobial Susceptibility Testing (EUCAST, 2000) to estimate the MICs of an antimicrobial agent. Based on these studies, the MIC for each strain was defined as the lowest concentration of CuCl₂ that consistently prevented visible growth in at least three independent assays. Based on this definition, our raw data were images of a solid medium with increasing millimolar (mM) concentrations of CuCl₂ with or without the growth of resistant or susceptible strains. Representative examples of growth monitored photographically at 2, 4, and 6 days post-propagation are shown in Fig. S2. In this example, the PPFM M. sp. able to grow in 2 mM CuCl₂ was excluded from this study because its genome sequence is not available. Numerous images were used to construct a table of permissive and inhibitory concentrations of CuCl₂ for each strain (Table S2). The workflow used to assess MICs is shown in Fig. S3.

Medium, CuCl₂ solutions, inoculation of plates, and reproducibility of MICs

To maintain the reproducibility of MICs in different assays we used agar plates with chemically defined medium instead of solid PY-rich medium. The variation observed in PY plates may be due to the high binding affinity of Cu(I) for amino acids.

To maintain repeatability among assays, the range of CuCl₂ concentrations evaluated (0–2 mM) was added to the mm from a 50 mM stock solution of CuCl₂-2H₂O (Sigma‒Aldrich, St Louis, MO, USA) prepared in Milli–Q water and filter sterilized. The stock solution of CuCl₂-2H₂O was used on the day of preparation and then discarded.

Preparation of the inoculum, in the same way, was also important to avoid variation among MIC assays. Each strain was propagated in solid PY medium. After 3 days at 30 °C, the cultured bacteria were inoculated into 3 ml of PY broth. After overnight growth, cultures were adjusted to OD₆₂₀ = 0.3, washed twice with 10 mM MgSO4, and tenfold serially diluted (10⁻¹–10⁻⁶). Twenty microliters from each dilution were spotted on solid Mm supplemented with increasing millimolar (mM) concentrations of CuCl₂.

Occurrence of copper translocating P_1B-type ATPases (Cu-ATPases) encoded in the genomes of PAFD and PPFM

Cu-ATPases encoded in the PAFD and PPFM genomes were searched by BLASTP (https://blast.ncbi.nlm.nih.gov/Blast.cgi). The sole Cu-ATPase encoded in the genome of R. etli CFN42 was used (WP_011427866.1) as a query protein, because of functional evidence previously demonstrated by the copper sensitivity phenotype of a mutant actP::ΩSp (González-Sánchez et al., 2018). Cu-ATPases containing eight transmembrane helices (TMHs) typical for Cu-ATPases were filtered from metal-ATPases using TOPCONS (http://topcons.net/). Additionally, we used CD-Search (https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi) to verify the presence of the P-type ATPase Cu-like domain cd02094 (NCBI) containing three invariant amino acid signature motifs involved in Cu(I) binding and translocation: two cysteine residues CXC in TMH6; one tyrosine, one asparagine, and one proline YN(X₄)P residue in TMH7; and one methionine followed by serine residue MXXSS in TMH8 (Arguello, 2003). The multiple sequence alignment of Cu-ATPases performed with Clustal Omega version 1.2.4 at MBL-EBI (https://www.ebi.ac.uk/Tools/msa/clustalo/) shows the signature motifs shared in 28 Cu-ATPases (Fig. S4). To verify that we did not miss Cu-ATPases by BLASTP searches, we retrieved all the metal translocating P-type ATPases annotated in the target genomes from the Bacterial and Viral Bioinformatics Resource Center (https://www.bv-brc.org/) formerly Pathosystems Resources Integration Center (https://patricbrc.org). The Cu-ATPases were filtered by searching the signature motifs with CD-Search. The occurrence analysis is summarized in Table 2. The 28 amino acid sequences obtained are shown in the Table S3 (Cu-ATPase sequences).

Evolutionary divergence among Cu-ATPases encoded in the genomes of PAFD and PPFM

The 28 amino acid sequences, obtained as described above, were used to infer their phylogeny. Four sequences were used as an outgroup: two putative zinc transporters encoded in the genomes of Mr. extorquens AM1 (WP_012752954.1) and Mb. organophilum DSM760 (PVZ05212.1), as well as two Cu-ATPases encoded in the genome of Acidithiobacillus ferrivorans ACH (WP_215894663.1) and (WP_215895207.1). Multiple sequence alignment, phylogenetic reconstruction analyses, and visualizations of phylogenetic data were performed using the pipeline of ETE3 at GenomeNet (https://www.genome.jp/tools/ete/). The workflow was as follows: clustalo_default-trimal001-prottest_default-phyml_default_boostrap. Branch support was computed using 100 bootstrapped trees. The Cu-ATPase sequences used to infer the phylogeny are available in Table S4 (Cu-ATPases for phylogeny).

Contextual information from conserved gene neighborhood (GN) analysis

This analysis was performed using TREND (Gumerov & Zhulin, 2021), a freely available bioinformatic tool for the tree-based exploration of neighborhoods and domains (http://trend.evobionet.com/). The five Cu-ATPases (CopA1-CopA5) encoded in the genome of Mr. extorquens AM1 (Table S5, CopA sequences for TREND) were used as input information in the neighborhoods pipeline. As a result, we obtained a phylogenetic tree combined with interactive gene neighborhoods that showed predicted operons based on the distance between genes, provided domains and features of gene products, and identified homologous proteins. Unexpectedly, the copA2. gene neighborhoods could not be obtained using the neighborhoods pipeline of TREND. These data were manually obtained at NCBI using CopA2 NCBI Refseq as input (WP_158022369.1) (see the workflow in Fig. S5).

Inference of the genome-based Cu efflux resistome in the PPFM Mr. extorquens AM1

The thoughtful curation of the contextual information obtained for the five copA genes was used to infer five different Cu-Homeostasis Gene Clusters (Cu-HGCs). Their gene organization is listed in Table S6, (Cu-HGCs). Their potential mobility was inferred by the presence of proteins whose predicted function has been associated with mobile genetic elements (MGEs) and horizontal gene transfer (HGT). The proteins encoded in each Cu-HGC are listed in Table S6, (Cu-HGCs).

Comparison of CuCl₂ MICs between PAFD and PPFM revealed hypertolerance of Mr. extorquens AM1

Based on the MIC of copper reported for bacteria isolated from mineral deposits, polluted environments and clinical samples we defined the range of high tolerance between 1- and 20-mM (Cusick et al., 2021; Watkin et al., 2009; Yik et al., 2018). In the present study, the MIC for each strain was obtained as described in the Material and Methods. Figure 1 shows the MIC of CuCl₂ that consistently prevented visible growth of the PAFD and PPFM strains in three independent assays. Figures 1A and 1B show that PPFM are more tolerant to CuCl₂ than PAFD. The highest CuCl₂ MIC observed in Mr. extorquens AM1 (1.9 mM) was similar to the MICs reported for bacterial strains isolated from mines, industrial wastewater, and clinical samples, such as S. meliloti CCNWSX0020 (1.8 mM), C. metallidurans CH34 (2 mM), and A. baumannii (1.5 mM) (Li et al., 2014; Wiesemann et al., 2013; Williams et al., 2016).

Several PAFD and PPFM harbor multiple Cu-ATPases

Copper extrusion from the cytoplasm to the periplasm is the most prevalent and crucial mechanism for maintaining the homeostasis of copper. It is performed mainly by P_1B-1-type ATPases known as CopA (Argüello, Raimunda & Padilla-Benavides, 2013). Since the presence of multiple CopA is a characteristic highly conserved in the genome of Rhizobiales (Cubillas et al., 2017), we performed an occurrence analysis to investigate the possibility of a correlation between the number of Cu-ATPases in the genomes of PAFD and PPFM and their MICs. Table 2 indicates that the number of Cu-ATPases per genome varied from 1 to 5, and their role in copper translocation was supported by the presence of amino acid motifs CXC, YN(X₄)P, and M(X₃)S located at transmembrane helices 6, 7, and 8, respectively (Fig. S4). Although the highest MIC of Mr. extorquens AM1 (1.9 mM) correlated with its highest number of Cu-ATPases (5), no clear correlation could be established between MICs and the number of Cu-ATPases in other species.

Cu-ATPases from PAFD and PPFM show high genetic diversity

Table 2 shows that multiple Cu-ATPases are encoded in different replicons, namely, chromosomes, megaplasmids, and plasmids. This result suggests that genes encoding Cu-ATPases may have different evolutionary histories. To test this hypothesis, the evolutionary divergence of the 28 Cu-ATPases (Table S4 Cu-ATPases for phylogeny) was inferred by maximum likelihood phylogenetic analysis. The 28 Cu-ATPases from PAFD and PPFM were sorted into eleven distantly related lineages (Fig. 2). Large divergence was also observed among the multiple Cu-ATPases encoded in the genome of several species, suggesting evolution by HGT. The phylogeny also showed that the Cu-ATPases MeAM3Ch and MeAM4Mp from Mr. extorquens AM1 located in lineage X were closely related, suggesting evolution by gene duplication.

Gene neighborhoods analysis of Cu-ATPases uncovers a complex and potentially mobile Cu efflux resistome in Mr. extorquens AM1

The high copper resistance exhibited by Mr. extorquens AM1 and the presence of five Cu-ATPases (CopA1-CopA5) encoded in its genome led us to define the Cu efflux resistome based on contextual information of copA1-copA5 genes as detailed above in Materials and Methods. This contextual information revealed that the five CopA belong to five different Cu-Homeostasis Gene Clusters (Cu-HGCs). Detailed schemes for each Cu-HGC, based on the interactive gene neighborhood diagrams displayed in TREND, are shown in Figs. 3–5.

The proteins encoded in each Cu-HGC are listed in Table S6 (Cu-HGCs). The organization of the Cu efflux resistome in five Cu-HGCs is described below.

Chromosomal 23.9 kb Cu/Ni/Co-homeostasis gene cluster (Cu/Ni/Co-HGC)

We first analyzed the GN of the Cu-ATPase coding gene copA3 due to the close phylogenetic relationship with CopA4, which suggests a gene duplication event (Fig. 2, lineage X). A careful inspection revealed that copA3 is part of a putative Cu-resistance operon (Fig. 3A, blue dashed line). Three additional metal resistance operons were predicted, with one upstream (light yellow dashed line) and two downstream of the Cu-resistance operon (Fig. 3A, dashed lines light yellow, and blue).

Cooccurrence of Cu-ATPases/CusAB and Cu-chaperones in the Cu-resistance operon

The first gene of this operon encodes Cu-ATPase CopA3, which is presumed to be a copper efflux pump (Fig. 3A, highlighted in yellow). Six hundred base pairs downstream of the copA3 gene is encoded a cytoplasmic short cysteine-rich protein with a length of 110 amino acids and a four-helical copper-binding bundle that might have a role in copper storage or as a Cu chaperone capable of binding Cu(I) and delivering it to the inner membrane CopA which transports Cu(I) to the periplasm. One thousand base pairs downstream of genes coding Cu-binding proteins are encoded two transporters belonging to the resistance-nodulation and cell division (RND) superfamily. The first shares 70% identity with the inner membrane permease CusA from E. coli, and the second shares 80% similarity with the periplasmic membrane fusion protein from E. coli CusB. In enterobacteria, both proteins are part of a cation antiporter CusCFBA system for Cu(I) and Ag(I) detoxification (Delmar, Su & Yu, 2015). However, neither the outer membrane factor CusC nor the periplasmic metallochaperone CusF homologs were predicted by GN analysis 10 kb upstream and downstream of cusAB genes. The polypeptide of 73 amino acids encoded downstream cusB gene (Fig. 3A) contains a heavy metal-associated domain shared among Cu(I) chaperones; however, two-sequence BlastP comparisons did not identify significant alignments with CusF. Moreover, PSORT (https://www.psort.org/) and BUSCA (http://busca.biocomp.unibo.it/), two web servers for subcellular localization of proteins, could not predict a periplasmic localization characteristic of E. coli CusF. Only the protein WP_003603728.1 produced significant alignment (25.9% identity, 97% query cover and 2e-²³ e-value) with CusC from E. coli (WP_000074194.1). Both shared domains with outer membrane factors operate in conjunction with membrane fusion proteins to transport substrates across membranes. This putative Mr. extorquens CusC is part of a multidrug efflux system together with one CusB and two CusA proteins located far from the Cu-resistance operon.

Putative Ni/Co resistance rcnA/rcnR operon

Downstream of the Cu-resistance operon, TREND predicted a couple of genes organized in a single operon with significant similarity (67%) to E. coli nickel/cobalt efflux proteins RcnA, and its transcriptional repressor, RcnR (Fig. 3A, light yellow dashed line).

Pb/Cd/Zn homeostasis operon

An operon harboring seven genes (Fig. 3A, dark blue dashed line) is located downstream of the rcnA/rcnR operon in the complementary DNA strand. Thorough curation of the function of these proteins did not reveal significant information on their role in Cu homeostasis.

Nucleases, CRISPR/Cas system, and mobility of Cu/Ni/Co-HGC

This multimetal homeostasis gene cluster is bordered at the 5′ end by a two-gene operon coding an exoDNAse V and a CRISPR/Cas system-associated protein Cas4 both sharing RecB domains (Fig. 3A, orange rectangles). The 3′ end is flanked by a type II endonuclease (orange rectangle), which is probably a remnant of the restriction-modification system (RM). It has been well documented that to maintain long-term persistence in the genome, RM and CRISPR/Cas systems are often linked with mobile genetic elements that confer an adaptive advantage. Their presence in the Cu/Ni/Co-HGC suggests that this gene cluster could be part of a genomic island prone to HGT (Oliveira, Touchon & Rocha, 2014).

Cu homeostasis gene cluster-1 (Cu-HGC-1, 25.7 kb)

Figure 2 suggests a close phylogenetic relationship between CopA3 and CopA4. The GN of CopA4 revealed that this megaplasmid located Cu-HGC-1 and the chromosomal Cu/Ni/Co-HGC described above in the section “Gene neighborhoods analysis of Cu-ATPases uncovers a complex and potentially mobile Cu efflux resistome in Mr. extorquens AM1” share not only their Cu-ATPases but also the Cu-resistance and cbb₃-Cox biogenesis operons (Fig. 3B, green and light brown dashed lines).

TREND also predicted an extra copy of CusA and a fragment of CusB that might be the result of a genomic rearrangement involved in the duplication of this region. Cu-HGC-1 is flanked at the 5′ and 3′ ends by transposases belonging to the IS3, IS110, and IS630 families that resemble the structure of composite transposons. In addition, a protein with an integrase domain is also located at the 5′ end of Cu-HGC-1. The presence of these enzymes suggests the potential intra- and intergenomic mobility of this gene cluster.

Megaplasmid-located Cu-homeostasis gene cluster-2 (Cu-HGC-2) encodes Cu efflux redundant functions

Figure 4A shows the gene neighborhood analysis of the Cu-ATPase CopA5. Genes located upstream and downstream of copA5 encode the RND transporters CusB and CusA, and the Cu chaperones CopZ and CusF, respectively. These genes span 16 kb and together integrate into the second Cu homeostasis gene cluster (Cu-HGC-2). No genetic elements involved in mobility were detected.

An unknown interrelation between efflux of copper and modulation of membrane charge highly conserved among the PPFM was found in the Cu-homeostatic gene cluster-3 (Cu-HGC-3)

The Cu-ATPase-encoding gene copA1 is located in a chromosomal locus. According to TREND, copA1 forms part of an operon together with two downstream genes (Fig. 4B). The product of the first gene shows significant identity (58%) with CueR from E. coli, a Cu(I)-responsive transcriptional regulator protein. The second gene encodes a putative M20 family metallopeptidase. Since almost all naturally occurring metallopeptidases are zinc-dependent enzymes, the functional relationship between the products of these three genes could not be established. Upstream of the copA1 gene, TREND predicted a second two-gene operon. The first gene, lpiA, encodes a polypeptide of 343 amino acids with a significantly similar identity (31.16%) with the lysylphosphatidylglycerol synthetase from R. tropici CIAT899 (WP_015341009.1) responsible for lysylphosphatidylglycerol (LPG) formation. LPG produces a positive charged membrane that defends bacteria from cationic antimicrobial molecules. This cytoplasmic membrane lipid is synthesized only at low pH (4.5) and is involved in acid tolerance (Sohlenkamp et al., 2007). Downstream lpiA, is encoded a sensor kinase member of the PhoR family (WP_012753105) that may be involved in sensing phosphate status. The functional interrelation between lpiA and copA is unknown, it can be speculated the efflux of copper under low pH conditions. This genetic organization is not shared between the acid-tolerant strain R. tropici CIAT899 and Mr. extorquens AM1 but is highly conserved in the PPFM studied in this work.

The mobile genetic element CopA2 (MGE-CopA2, 24.7 kb) shares characteristics of both integrative and conjugative elements (ICE), and composite transposons (CTn)

The gene encoding the Cu-ATPase CopA2, along with numerous genes often present in MGEs, are located in a chromosomal locus (Fig. 5). In the complementary DNA strand, upstream copA2 gene is encoded a site-specific integrase that may catalyze the integration/excision of MGE. Downstream of the integrase gene, is encoded a type II toxin-antitoxin system (highlighted in blue), which is important for the long-term persistence of MGEs in their hosts (van Melderen, de Bast & Rosenberg, 2009). Downstream this T-A system is encoded an antirestriction protein homologous to ArdC (highlighted in gray), which may inhibit type I and type II R-M systems of the recipient cell and avoid degradation of the incoming MGE upon an HGT event (González-Montes et al., 2020). The presence of two IS3 family transposases bordering this genetic region (orange rectangles) resembles the structure of a composite transposon. The complexity of this putative MGE increases due to the presence of a traG gene at the 5′ end (highlighted in yellow). The product of this gene may be a component of the type IV secretory system essential for DNA transfer in bacterial conjugation. The relaxase domain-containing protein (highlighted in yellow) encoded four kb downstream of traG may initiate the conjugative transfer of DNA binding to the origin of transfer (oriT) and melt the double helix. The copA2 gene, together with the integrative and conjugative transfer-associated functions described above, resembles the structure of a complex MGE known as an integrative and conjugative element (ICE). An ICE exhibits two different states: an integrative state, in which its DNA resides in the chromosome of the host, and a conjugative state, in which its DNA is excised from the chromosome of the host and can potentially spread horizontally by conjugative transfer to a new cell (Delavat et al., 2017).