Genomic identification and characterization of the elite strains Bradyrhizobium yuanmingense BR 3267 and Bradyrhizobium pachyrhizi BR 3262 recommended for cowpea inoculation in Brazil

Abstract

The leguminous inoculation with nodule-inducing bacteria that perform biological nitrogen fixation is a good example of an “eco-friendly agricultural practice”. Bradyrhizobium strains BR 3267 and BR 3262 are recommended for cowpea (Vigna unguiculata) inoculation in Brazil and showed remarkable responses; nevertheless neither strain was characterized at species level, which is our goal in the present work using a polyphasic approach. The strains presented the typical phenotype of Bradyrhizobium with a slow growth and a white colony on yeast extract-mannitol medium. Strain BR 3267 was more versatile in its use of carbon sources compared to BR 3262. The fatty acid composition of BR 3267 was similar to the type strain of Bradyrhizobium yuanmingense; while BR 3262 was similar to Bradyrhizobium elkanii and Bradyrhizobium pachyrhizi. Phylogenetic analyses based on 16S rRNA and three housekeeping genes placed both strains within the genus Bradyrhizobium: strain BR 3267 was closest to B. yuanmingense and BR 3262 to B. pachyrhizi. Genome average nucleotide identity and DNA–DNA reassociation confirmed the genomic identification of B. yuanmingense BR 3267 and B. pachyrhizi BR 3262. The nodC and nifH gene analyses showed that strains BR 3267 and BR 3262 hold divergent symbiotic genes. In summary, the results indicate that cowpea can establish effective symbiosis with divergent bradyrhizobia isolated from Brazilian soils.

Introduction

Bacteria collectively known as rhizobia form an important part of the soil microbiota and perform biological nitrogen fixation (BNF) through nitrogenase activity when in symbiosis with leguminous plants. This ecological phenomenon has great biotechnological impact on biomass and grain production. For example, in the soybean crop production in Brazil it is estimated to save over US$ 10 billion annually by the use of Bradyrhizobium inoculation instead of chemical fertilization.¹

The Brazilian Ministry of Agriculture has a list of rhizobial strains recommended for more than 50 leguminous grain-producing crops, forage and green manure.² In recent years efforts have been made to better characterize these strains recommended for inoculation in Brazil and at least five new species have been described.^3–6

Brazil is the world's third leading cowpea producer, with an estimated production of 500,000 tons per year.⁷ The crop yield varies from 400 to 2000 kg ha⁻¹, depending on the system and the region of cultivation.⁷ This yield has been rising in recent years with improvements in cultivation management, such as the inoculation of seeds with nitrogen-fixing bacteria, a practice applied to approximately 100,000 hectares.

The strains BR 3267 and BR 3262 are considered as “elite” for inoculation of cowpea plants in Brazil² and several studies under controlled and field conditions have shown that both strains make significant contributions to crop yields, including more than 50% N accumulation via BNF.⁸ BR 3267 strain was isolated from the semiarid northeastern region of the country, using cowpea as trap plants, while BR 3262 strain was isolated from an Atlantic Forest area in southeastern Brazil using the same strategy.^1,9 Although these strains were isolated more than a decade ago, they were only partially characterized through assessment of the growth rate, colony morphology and 16S rRNA phylogeny,¹¹ but not classified at the species level.

Zilli and colleagues¹¹ have characterized the 16S rRNA genes of the BR 3267 and BR 3262 strains and concluded that both are members of the genus Bradyrhizobium. This genus was created in the early 1980s to accommodate root nodule-inducing bacteria with slow growth on media containing mannitol and yeast extract.¹² Since the turn of the century, two major subgroup divisions (I and II) have been recognized within this genus based on DNA–DNA hybridization.^12–14 The Bradyrhizobium japonicum and Bradyrhizobium elkanii were the first species assigned to subgroup I and II, respectively.^15,16 According to Zilli and colleagues,¹⁰ BR 3267 clustered within the subgroup B. japonicum and BR 3262 within the subgroup B. elkanii. However, in the light of current knowledge, those results can be considered inconclusive because new Bradyrhizobium species was recently described. Thus, further investigation employing the latest molecular techniques is needed for the correct positioning of these strains.

In the past five years, the use of housekeeping genes as powerful phylogenetic markers for bacteria has led to the description of over 15 new species within the genus Bradyrhizobium¹⁷ and enabled the separation of genetically close strains into different species. Examples are the definition of the species Bradyrhizobium diazoefficiens based on B. japonicum and Bradyrhizobium pachyrhizi from B. elkanii.^3,18 Furthermore, new methods for genome-to genome comparison, like Average Nucleotide Identity (ANI) and Genome Blast Distance Phylogeny (GBDP) have been introduced and are contributing to improve bacterial taxonomy.^19,20

Therefore, the goal of our study was to characterize at a finer the taxonomic level both BR 3267 and BR 3262 strains using a polyphasic approach.

Materials and methods

Strains used

The cowpea strains BR 3267 and BR 3262, and the type strain B. elkanii USDA 76^T were obtained from the Johanna Döbereiner Biological Resource Center (CRB-JD, Embrapa Agrobiologia, Seropédica-Rio de Janeiro, Brazil). The type strains B. pachyrhizi PAC 48^T (=LMG 24246^T) and Bradyrhizobium yuanmingense CCBAU 10071^T (=LMG 21827) were obtained from the LMG Culture Collection (Belgium). The strains were grown on yeast extract-mannitol agar medium (YMA) and were incubated at 28 °C²¹ for seven days until they reached sufficient colony growth levels to observe morphological features and purity.

Phenotypic and physiologic characterizations

Inoculum preparation for both BR 3267 and BR 3262 strains was carried out using YMA medium at 28 °C. Carbon source utilization was assessed with Biolog GN2 microplates (Biolog Inc., Hayward, CA) following the manufacturer's instructions, except that cell concentration was adjusted to 5 on the McFarland scale. Plates were incubated in the dark at 28 °C for ten days. The biochemical features were assessed using API 20 NE strips (bioMérieux, Marcy-L’Etoile, France) following a standard protocol that uses a saline solution (0.85% NaCl) for bacterial suspension. Additionally, the tolerance to abiotic stress, such as temperature, pH and salinity (NaCl), was determined by examining the growth in YMA medium. The temperature tolerance was evaluated at 15, 20, 25, 28, 30, 32 and 37 °C, and the pH tolerance was tested in a range from 4 to 10. The salinity tolerance was examined at 28 °C in YMA medium supplemented with 0.1, 0.3, 0.5, 1.0, 1.5, 2.0 or 2.5% (w/v) of NaCl. The resistance to antibiotics was determined using YMA medium and the disk diffusion method for ampicillin (25 μg), chloramphenicol (50 μg), erythromycin (30 μg), gentamicin (10 μg), kanamycin (30 μg), neomycin (10 μg), penicillin (10 μg), streptomycin (10 μg) and tetracycline (30 μg). All tests were run in triplicate.

Fatty acid composition

The BR 3267 and BR 3262 and type strains B. elkanii USDA 76^T, B. pachyrhizi PAC 48^T and B. yuanmingense CCBAU 10071^T type strains were characterized based on whole-cell fatty acids, derived using the methyl esters’ method (FAME) and analyzed with a Hewlett Packard gas chromatograph fitted with a fused silica capillary column (25 m × 0.2 mm internal diameter). The type strains B. elkanii USDA 76^T, B. pachyrhizi PAC 48^T and B. yuanmingense CCBAU 10071^T were chosen because they were the closest related strains to BR 3262 and BR 3267 respectively, based on EzTaxon identification search (http://www.ezbiocloud.net/eztaxon/identify). ChemStation A.09.01 software [1206] and a MIDI Microbial Identification System 4.0 (Sherlock TSBA Library, MIDI ID, Inc., Newark, DE, USA) were used for the fatty acid identification.

DNA extraction, PCR amplification and sequencing of PCR fragments

Pure cultures of Bradyrhizobium strains were grown on YM medium for 4 days at 28 °C under stirring of 120 rpm. Then, 2 mL of cell suspension were centrifuged (12,000 × g), and the total genomic DNA was extracted using a Wizard^® Genomic DNA Purification Kit (Promega, USA) following the manufacturer's instructions. The recA gene was amplified by PCR using the primers TSrecAf and TSrecAr, and the glnII gene was amplified with the primers TSglnIIf and TSglnIIr.²² The primers NodCfor540 and NodCrev1160 were used to amplify the symbiotic nodC gene, as described by Sarita et al.²³ For the nifH, the primers NifHF/NifHI were used following the methods described by Laguerre et al.²⁴ DNA sequences of the PCR products were obtained for both strands with the same primers used for gene amplifications. The 16S rRNA sequences of the BR 3267 (Gene Bank A.no. AY649439) and BR 3262 (A.no. AY649430) strains obtained in a previous study¹⁰ were retrieved from GenBank (http://www.ncbi.nlm.nih.gov/genbank/). For the gyrB gene, the sequences were retrieved from draft genome sequences of the BR 3267 and BR 3262 strains (Ac. no. KT005410 and KT005409, respectively).

Phylogenetic analysis

Forward and reverse readings of the 16S rRNA, recA, glnII, nodC and nifH gene fragments were edited and assembled using DNABASER (http://www.dnabaser.com). Multiple sequence alignments were performed using ClustalW through MEGA 6.0.²⁵ Maximum likelihood (ML) phylogenetic reconstructions were completed using MEGA 6 for single gene sequences (16S rRNA, recA, glnII, gyrB, nodC and nifH), and a multilocus sequence analysis (MLSA) was performed using the concatenated sequences of the housekeeping genes recA, glnII and gyrB. The distance matrices were calculated using the Kimura two-parameter substitution model,²⁶ and the robustness of the tree nodes was evaluated with a bootstrap analysis²⁷ using 500 pseudoreplicates. The software default parameters were considered in all of the analyses. Sequences of rhizobial type strains used for alignment and phylogenetic analysis were retrieved from GenBank (http://www.ncbi.nlm.nih.gov/genbank/), and a Microvirga vignae BR 3299 strain²⁸ was used as the outgroup in the ML analysis of the 16S rRNA gene.

DNA–DNA relatedness and average nucleotide identity (ANI)

The DNA–DNA hybridization (DDH) among bacterial strains was determined based on the thermal denaturation temperatures of hybrid and homologous genomic DNA, as described by Gonzalez and Saiz-Jimenez,^29,30 except that the DNA concentration was adjusted to 2 μg per reaction. The experiments were performed in 96-well optical plates with appropriate optical adhesives in three replicates and including wells without DNA as negative control. For DNA hybridization, a Bio-Rad MyCicler thermocycler was used, and for fluorimetric measurements during the denaturation stage an Applied Biosystems 7500 Real-Time system was used. Thermal conditions for hybridization consisted on a denaturation step of 99 °C for 10 min, followed by an annealing period of 8 h at 79° C (optimum temperature for renaturation – TOR).^29,31 It was followed by progressive 10-min steps, each at 1.8 °C below the previous one, until 25 °C when it was hold for 30 min before refrigeration to 4 °C. The fluorescence was measured using a denaturation ramp settled at the step and hold mode. Heating rate was 0.2 °C s⁻¹ with fluorescence decreasing measurement at each 0.2 °C step, during a 12 s hold, and between 25 and 99.9 °C.²⁹ The DDH experiments were performed between strains BR 3267 and BR 3262 and their phylogenetically related type strains based on the 16S rRNA gene similarities. BR3267 was compared against B. yuanmingense CCBAU 10071^T, and BR 3262 was compared against B. pachyrhizi PAC 48^T and B. elkanii USDA 76^T.

The ANI estimation was performed with a JSpecies platform version 1.2.1,³² and MUMmer was used for genome alignment.³³ All system requirements (BLAST and MUMmer) were downloaded and installed locally in a system Linux. Bacterial genome were retrieved from NCBI database as fasta file and used to ANI calculation, considering JSpecies default parameters. The sequenced genome of strain BR 3262 (Ac.no. LJYE00000000.1)³⁴ was compared to genomes of the strains USDA 76^T (B. elkanii – GenBank A.no. ARAG00000000.1) and PAC48^T (B. pachyrhizi – GenBank A.no. LFIQ00000000.1). The genome of the strain BR 3267 (A.no. LJYF00000000.1)³⁵ was compared to that of the strain CCBAU 35157 (B. yuanmingense – A.no. AJQL00000000.1).

Results

Phenotypic and physiologic features

Both BR 3267 and BR 3262 strains can grow in a pH range of 4–10, tolerate up to 0.5% NaCl in the medium and present normal growth in a temperature range from 15 to 32 °C, and BR 3267 could also grow up to 37 °C (Table 1). Likewise, both strains showed positive reactions to enzyme urease and nitrate reductase. Their susceptibilities to different antibiotics were variable (Table 1). BR 3267 was sensitive to streptomycin and tetracycline, while BR 3262 was tolerant. With respect to the carbon sources evaluated with the Biolog kit, BR 3262 strain was able to use 33 of the 95 sources tested, while BR 3267 strain was able to use more than 50.

Table 1.

Phenotypic features of the cowpea Bradyrhizobium strains BR 3267 and BR 3262 as characterized by API 20NE and Biolog GN2 microplates. Data obtained from 5 days duplicate read mean values (+ = positive; − = negative, w = weak).

Phenotypic feature	BR 3267	BR 3262	Phenotypic feature	BR 3267	BR 3262	Phenotypic feature	BR 3267	BR 3262
Enzymatic reaction			pH 4	+	+	Ester
Catalase	+	−	pH 10	+	+	Methyl pyruvate	+	−
Hydrolysis of gelatine	+	−	1% NaCl	−	−	Mono-methyl-succinnate	+	−
Hydrolysis of esculin	+	+	Carbohydrates			Amino acids
β-galactosidase	w	−	d-Arabitol	+	−	d-Alanine	−	+
Resistance to (μg)			α-d-Glucose	+	−	l-Alanine	−	+
Ampicillin (25)	−	+	d-Mannose	+	−	Glycyl-l-aspartic acid	+	−
Chloramphenicol (50)	−	+	d-Psicose	+	−	l-Leucine	+	−
Kanamycin (30)	+	−	d-Sorbitol	+	−	l-Phenylalanine	+	−
Neomycin (10)	+	−	Carboxylic acids			l-Proline	+	−
Penicillin (10)	−	+	d-Galacturonic acid	+	−	d-Serine	+	−
Streptomycin (10)	−	+	d-Glucosaminic acid	+	−	l-Threonine	+	−
Tetracycline (30)	−	+	α-Ketovaleric acid	+	−	γ-Aminobutyric acid	+	−
Erythromycin (30)	+	+	d,l-Lactic acid	+	−	Polymer
Gentamicin (10)	+	+	Malonic acid	+	−	Dextrin	−	+
Growth at			Quinic acid	+	−	Amide
15 °C	+	+	d-Saccharic acid	+	−	Glucuronamide	+	−
32 °C	+	+	Sebacic acid	+	−	Amine
37 °C	+	−	Succinic acid	+	−	Phenylethylamine	+	−

Both strains BR 3267 and BR 3262 were positive to: Oxidase, Urease, Nitrate reduction, Hydrolysis of esculin, l-Arabinose, l-Fucose, d-Galactose, d-Mannitol, l-Rhamnose, Tween 40, Tween 80, Acetic acid, Citric acid, Formic acid, d-Galactonic acid lactone, d-Gluconic acid, α-Hydroxybutyric acid, β-Hydroxybutyric acid, γ-Hydroxybutyric acid, p-Hydroxyphenylacetic acid, α-Ketobutyric acid, α-Ketoglutaric acid, Propionic acid, Glycerol, l-Alanyl-glycine, l-Asparagine, l-Aspartic acid, l-Glutamic acid, Glycyl-l-glutamic acid, l-Pyroglutamic acid, l-Serine, Succinamic acid, l-Alaninamide and Bromo succinic acid; and negative to: Tryptophan deaminase, Glucose fermentation, Arginine dihydrolase, N-Acetyl-d-galactosamine, N-Acetyl-d-glucosamine, Adonitol, d-Cellobiose, i-Erythritol, d-Fructose, Gentiobiose, m-Inositol, α-d-Lactose, Lactulose, Maltose, d-Melibiose, β-Methyl-d-Glucoside, d-Raffinose, Sucrose, d-Trehalose, Turanose, Xylitol, α-Cyclodextrin, Glycogen, cis-Aconitic acid, d-Glucuronic acid, Itaconic acid, 2,3-Butanediol, d,l-α-Glycerol phosphate, Glucose-1-phosphate, Glucose-6-phosphate, l-Histidine, Hydroxy-l-proline, l-Ornithine, Urocanic acid, Inosine, Uridine, Thymidine, Putrescine and 2-Aminoethanol.

Fatty acid composition

The fatty acid composition was analyzed for the BR 3262 and BR 3267 strains and for B. elkanii USDA 76^T, B. pachyrhizi PAC 48^T and B. yuanmingense CCBAU 10071^T, the closest type strains based on EzTaxon Identify search of the 16S rRNA gene (http://www.ezbiocloud.net/eztaxon/identify). The fatty acid compositions of BR 3267 and BR 3262 were similar to those of other strains within Bradyhizobium genus, especially because of the presence of 16:00 in approximately 15% and a summed feature 8 (18:1 w7c) between approximately 70 and 80%.⁴⁰ BR 3262 has a composition close to those of B. elkanii USDA 76^T and B. pachyrhizi PAC 48^T, with the presence of 19:0 cyclo w8c of 15.98% and a summed feature 8 (18:1 w7c) of 69.32. Instead, BR 3267 and B. yuanmingense CCBAU 10071^T did not present 19:0 cyclo w8c and had a summed feature 8 of 86.42 and 84.63%, respectively.

Phylogenetic analyses

The Maximum Likelihood (ML) phylogenetic analysis performed on the 16S rRNA (1227 bp) gene placed the cowpea strains BR 3267 and BR 3262 into subgroups I and II, respectively, within the genus Bradyrhizobium (Fig. 1). The closest strain to BR 3267 was B. yuanmingense CCBAU 10071^T (99.5% 16S rRNA similarity), with 67% ML bootstrap support, followed by Bradyrhizobium subterraneum 58 2-1^T (99.6% 16S rRNA similarity) in a separated branch with 62% bootstrap. In contrast, BR 3262 was placed within subgroup II of Bradyrhizobium in a phylogenetic linage together with B. pachyrhizi PAC48^T (100% 16S rRNA similarity), B. elkanii USDA 76^T (99.9%) and Bradyrhizobium tropiciagri CNPSo 1112^T (99.2%), but with nodes lower than <60% of bootstrap. The Bradyrhizobium ferriligni CCBAU 51502^T was the closest taxon to the clade of B. tropiciagri–B. elkanii–B. pachyrhizi – BR 3262, as shown by the high ML bootstrap support (91%, Fig. 1); it shares 97.8% of 16S rRNA similarity with BR 3262.

Fig. 1.

Maximum Likelihood phylogenetic tree based on 16S rRNA gene sequences showing relationships between the strains BR 3262 and BR 3267 and the type species (T) of the genus Bradyrhizobium. Bootstrap values were inferred from 500 replicates and are indicated at the tree nodes when ≥50%. GenBank accession numbers are provided in parenthesis. Microvirga vignae BR3299^T was included as the outgroup. The bar represents two estimated substitutions per 100 nucleotide positions.

We used aligned sequences from the recA (359 bp), glnII (505 bp) and gyrB (552 bp) housekeeping genes, which resulted in a concatenated sequence of 1416 bp, to better assess the phylogenetic affiliation of the BR 3267 and BR 3262 strains. ML phylogenetic trees of single and concatenated genes (Fig. 2) were generated and showed high congruence; therefore, we presented only the concatenated one. The concatenated phylogenetic analysis revealed a clear affiliation between BR 3267 and B. yuanmingense CCBAU 10071^T, with high confidence bootstrap values of the 100%, in a lineage that was distinct from the closest soybean isolates Bradyrhizobium daqingense CCBAU 15774^T48 and Bradyrhizobium huanghuaihaiense CCBAU 23303^T.⁴⁹ On the other hand, BR 3262 was phylogenetically close to B. pachyrhizi PAC48^T, which formed a lineage with high confidence values (bootstrap ML of 98%) and separated from the closed taxa B. ferriligni CCBAU 51502^T and B. elkanii USDA76^T. Additionally, pairwise comparisons of the concatenated sequences revealed that BR 3267 shared 99% sequence identity with B. yuanmingense CCBAU 10071^T, followed by Bradyrhizobium liaoningense USDA 3622^T (95%) and B. daqingense CCBAU 15774^T (94.5%). Likewise, BR 3262 was similar to B. pachyrhizi PAC48^T (98.6%), followed by B. elkanii USDA 76^T (96.9%) and B. ferriligni CCBAU 51502^T (95.5%). In an independent analysis of the recA, glnII and a concatenated recA-glnII sequences (data not shown), BR 3267 and BR 3262 were far from the newest proposal species Bradyrhizobium kavangense 14-3^T,⁵⁰Bradyrhizobium vignae 7-2^T51 and B. subterraneum 58-2-1^T52 isolated from African soils.

Fig. 2.

Unrooted maximum likelihood phylogenetic tree based on three concatenated sequences (recA, glnII, gyrB) showing relationships between the strains BR 3262 and BR 3267 and type species (T) of the genus Bradyrhizobium. The phylogenetic tree was built using the Kimura two-parameter method. Bootstrap values were inferred from 500 replicates and are indicated at the tree nodes when ≥50%. GenBank accession numbers are provided in parenthesis. The bar represents two estimated substitutions per 100 nucleotide positions.

DNA–DNA relatedness and average nucleotide identity

After obtaining the phylogenetic results, we established the contrasts for the analysis of the DNA–DNA homology following the protocol described by Gonzales and Saiz-Jimenez.²⁹ BR 3262 was compared with B. elkanii USDA 76^T and B. pachyrhizi PAC 48^T (the closest strains according to MLSA), and BR 3267 was compared with B. yuanmingense CCBAU 10071^T (also the closest strain in the MLSA). The hybridization of BR 3262 × B. elkanii USDA 76^T and B. elkanii USDA 76^T × B. pachyrhizi PAC 48^T presented ΔT_m values equal to 10.1 °C and 3.9 °C, respectively, and the hybridization of BR 3267 × B. yuanmingense CCBAU 10071^T presented a value of 3.9 °C (Fig. 3).

Fig. 3.

Melting curves generated with Applied Biosystems 7500 Real-Time System for ΔT_m determination between homologous and hybrid DNA. (a) Comparison between strain BR 3267 and B. yuanmingense CCBAU 1007^T; (b) between BR 3262 and B. pachyrhizi PAC 48^T; (C) between BR 3262 and B. elkanii USDA 76^T.

For the ANI calculation, we used the available genomes of B. elkanii USDA 76^T and B. pachyrhizi PAC 48^T for comparison with BR 3262 and, B. yuanmingense CCBAU 35157⁵⁴ for comparison with BR 3267 because no genome was available for the type strain of B. yuanmingense CCBAU 10071^T. The ANI between strain BR 3262 and B. elkanii USDA 76^T and between strain BR 3262 and B. pachyrhizi PAC 48^T was 94.6 and 95.3, respectively. On another hand the ANI in the comparison of BR 3267 and B. yuanmingense CCBAU 35157 was 96.5.

Phylogeny based on the nodC and nifH symbiotic genes

Maximum Likelihood phylogenetic analyses of the symbiotic genes nodC and nifH separated the BR 3267 and BR 3262 strains into two divergent phylogenetic lineages (Fig. 4A for nodC, Fig. 4B for nifH). For both nodC and nifH, BR 3267 formed an independent branch that clustered together with B. yuanmingense CCBAU 10071^T and the soybean type strains that represent the symbiovar glycinearum. In contrast, BR 3262 was grouped in an early branch linked to the related strains B. pachyrhizi PAC48^T, B. elkanii USDA 76^T and B. ferriligni CCBAU51502^T for both symbiotic genes (Fig. 4A and B).

Fig. 4.

Unrooted maximum likelihood phylogenetic tree based on nodC (A) and nifH (B) genes showing relationships between the strains BR 3262 and BR 3267 and type species (T) and reference strains of the genus Bradyrhizobium. The phylogenetic tree was built using the Kimura two-parameter method. Bootstrap values were inferred from 500 replicates and are indicated at the tree nodes when ≥50%. GenBank accession numbers are provided in the parenthesis. The bar represents five or two estimated substitutions per 100 nucleotide positions.

Discussion

The analyses in the culture medium showed that both strains, BR 3267 and BR 3262, have similar behaviors in comparison to other members of the Bradyrhizobium genus regarding tolerance to temperature, pH and NaCl concentration in the culture medium. These traits indicate they are tolerant to different edapho-climatic factors and make them well adapted to the conditions of the tropical soils from which they were isolated. This adaptability is corroborated by their good performance in cowpea inoculations.^8,36,37 The BR 3267 and BR 3262 strains were positive to urease and nitrate reductase enzymes reactions, which are apparently common characteristics among strains of the genus Bradyrhizobium.^38,39 It seems to be a common characteristic among Bradyrhizobium species of subgroup I (close to B. japonicum) to be sensitive to streptomycin and tetracycline, as showed by BR 3267. In contrast, BR 3262 was tolerant to these antibiotics, a common characteristic of Bradyrhizobium species belonging to subgroup II (close to B. elkanii).⁴⁰ The 50 carbon sources used by BR 3267 points to greater access to different soil carbon compounds. The results of the fatty acid composition analysis corroborate previous findings that indicate that the main differences between subgroups I and II of Bradyrhizobium are the 19:0 cyclo w8c and the summed feature 8 (18:1 w7c), which indicate that BR 3267 belongs to the subgroup of B. japonicum and BR 3262 to the subgroup of B. elkanii.^40,41

When the 16S rRNA genes of strains BR 3267 and BR 3262 were first studied,¹¹ it was observed that although BR 3267 could be grouped with B. japonicum, there were differences with respect to the B. japonicum type strain pattern. At the time, however, the species B. yuanmingense, which also belongs to subgroup I of Bradyrhizobium, had not yet been described; thus, there was no species indication. It was observed on the other hand that BR 3262 had a strong resemblance to B. elkanii by the 16S rRNA gene, which supported that the strain should be assigned to that species. At that time, the species B. pachyrhizi, which in terms of base composition of the 16S rRNA gene differs from strains B. elkanii USDA 76^T only on the order of 0.1%, had not yet been described. 16S rRNA sequence analysis is not always fine enough for species assignment^42,43; in the genus Bradyrhizobium some species share almost identical 16S rRNA gene sequences. Examples are the type strains of the species B. japonicum and B. liaoningense,³ along with B. elkanii and B. pachyrhizi.¹⁸ Due to the high conservation of the 16S rRNA gene among members of the Bradyrhizobium genus,^14,15 several authors have recommended an MLSA of protein-encoding (i.e., housekeeping) genes.^6,22,44–47 In our case, the MLSA results were congruent to the 16S rRNA analysis and clearly indicate that BR 3267 and BR 3262 as members of the species B. yuanmingense and B. pachyrhizi, respectively (Fig. 2). The DNA–DNA homology results confirm the results obtained in the analysis of housekeeping genes, i.e., BR 3262 belongs to the B. pachyrhizi species, and BR 3267 is a member of the B. yuanmingense species. As previously shown, a ΔT_m value equal to 5 °C can be used as the threshold to delineate a species, where lower values indicate equal species and higher values indicate different species.²⁹ A ΔT_m value of 5 °C corresponds to a homology value of approximately 70% obtained by traditional DNA–DNA hybridization techniques.^29,30 Moreover, we calculated the ANI, a robust method that replaces the traditional DNA–DNA hybridization techniques.¹⁹ This method involves comparing the sequences of two bacterial genomes, and a value of 95% is considered as the threshold.^19,53 In this case, when the homology is greater than 95%, the two strains belong to the same species; if the homology is lower, they belong to different species. Thus, the results of the comparisons again confirmed that BR 3262 belongs to B. pachyrhizi (ANI greater than 95%) and the same BR 3267 to B. yuanmingense (ANI greater than 95%).

Symbiotic genes related to nodulation (nod) and nitrogen fixation (nif) processes hold no taxonomic information for species definition.¹⁷ However, in our study, they showed good phylogenetic correspondence to the core genome ones. Interestingly, nodC and nifH gene phylogenetic analyses have been used to predict the host range nodulation capacity and for the delineation of symbiovars, i.e., groups of rhizobial strains with similar symbiotic behaviors concerning nodulation and nitrogen fixation capacity with legumes.⁵⁵ For the nodC gene, for example, a BLAST search in the NCBI database (http://blast.ncbi.nlm.nih.gov/) showed the Arachis hypogae strains CCBAU 53380 (Ac. no. KC509232) and CCBAU 51663 (Ac. no. KC509217), isolated in China, as the Bradyrhizobium strains that had the closest similarity (98%) with BR 3267, indicating that peanut might be a possible host for this latter strain. Instead, in the nodC gene, BR 3262 had 95% similarity with the cowpea strain Bradyrhizobium sp. VUMPE10 (Fig. 4A), isolated from an European soil and proposed as a new symbiovar (sv. vignae) within the Bradyrhizobium genus.⁵⁶ These data indicate that cowpea nodulates with very effective Bradyrhizobium strains that are holding divergent symbiotic genes which might represent distinct symbiovars.

Cowpea can nodulate with a wide range of rhizobia genotypes, Bradyrhizobium being the major symbiont group. Thies and colleagues⁵⁷ proposed that cowpea's promiscuity for nodulation with the bradyrhizobial population may have a limit. However, a global Bradyrhizobium collection for cowpea nodulation has not yet been well analyzed to identify true cowpea symbionts in light of the symbiovar concept. It is not well understood whether any Bradyrhizobium group is specialized to produce nodules and to fix nitrogen specifically in cowpea.

Conclusion

The results of the identifications of B. yuanmingense BR 3267 and B. pachyrhizi BR 3262 can be considered reliable. Furthermore, analysis of the nodC and nifH genes indicates that these two strains are similar to strains isolated on other continents. Strain BR 3262 possibly belongs to the vignae symbiovar, which was recently assigned for strains isolated from cowpea in Europe. Both strains hold similar symbiotic genomes with Bradyrhizobium strains outside of South America.

Apparently, cowpea is able to establish effective symbiosis with divergent bradyrhizobia strains isolated from Brazilian soils, confirming its promiscuous capacity for nodulation.

Conflicts of interest

The authors declare they have no conflicts of interest.

Acknowledgements

This study was supported by the Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq). We would like to thank the Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (Capes) for the scholarships, likewise to CNPq for the fellowship of research productivity (PQ).