ORIGINAL RESEARCH article

Front. Microbiol., 11 November 2021

Sec. Microbial Symbioses

Volume 12 - 2021 | https://doi.org/10.3389/fmicb.2021.754837

Identification of the Important Genes of Bradyrhizobium diazoefficiens 113-2 Involved in Soybean Nodule Development and Senescence

  • 1. Key Laboratory of Biology and Genetic Improvement of Oil Crops, Ministry of Agriculture and Rural Affairs of PRC, Oil Crops Research Institute of Chinese Academy of Agriculture Sciences, Wuhan, China

  • 2. Shenzhen Branch, Guangdong Laboratory of Lingnan Modern Agriculture, Genome Analysis Laboratory of the Ministry of Agriculture and Rural Affairs, Agricultural Genomics Institute at Shenzhen, Chinese Academy of Agricultural Sciences, Shenzhen, China

Article metrics

View details

8

Citations

3,5k

Views

1,6k

Downloads

Abstract

Legume nodule development and senescence directly affect nitrogen fixation efficiency and involve a programmed series of molecular events. These molecular events are carried out synchronously by legumes and rhizobia. The characteristics and molecular mechanisms of nitrogen fixation at soybean important developmental stages play critical roles in soybean cultivation and fertilizer application. Although the gene expression of soybean were analyzed in nodules at five important soybean developmental stages, information on the expression of rhizobial genes in these nodule samples is limited. In the present study, we investigated the expression of Bradyrhizobium diazoefficiens 113-2 genes in the nodule samples from five developmental stages of soybean (Branching stage, flowering stage, fruiting stage, pod stage and harvest stage). Similar gene expression patterns of B. diazoefficiens 113-2 were existed during optimal symbiotic functioning, while different expression patterns were found among early nodule development, nitrogen fixation progress and nodule senescence. Besides, we identified 164 important different expression genes (DEGs) associated with nodule development and senescence. These DEGs included those encoding nod, nif, fix proteins and T3SS secretion system-related proteins, as well as proteins involved in nitrogen metabolism, ABC transporters and two-component system pathways. Gene Ontology, KEGG pathway and homology analysis of the identified DEGs revealed that most of these DEGs are uncharacterized genes associated with nodule development and senescence, and they are not core genes among the rhizobia genomes. Our results provide new clues for the understanding of the genetic determinants of soil rhizobia in nodule development and senescence, and supply theoretical basis for the creation of high efficiency soybean cultivation technology.

Introduction

The symbiotic relationship between legume and soil rhizobia leads to the formation of root nodule and fixation of atmospheric nitrogen (Ferguson et al., 2010). Depending on the persistence of the nodule meristem and morphology characteristics of nodules, nodules in legumes can be divided into determinate nodules and indeterminate nodules (Hirsch, 1992). Nodules formed in Glycine max and Lotus japonicus are determinate nodules, they have non-persistent meristems and predominately contain one specific form of the bacterial development, which is determined by the nodule’s age (Karr and Emerich, 1988; Karr et al., 1990). Indeterminate nodules, such as those in Medicago truncatula and Astragalus sinicus, contain persistent meristems and a mixture of bacteroid forms in all stages of development (Vasse et al., 1990). In both determinate and indeterminate nodules, nodule development directly affects nitrogen fixation activity, and the premature senescence of nodules negatively regulates the efficiency of symbiotic nitrogen fixation (Serova et al., 2018). Delaying nodule senescence is a practical measure to increase nitrogen fixation amount and solve the problem of late defertilization in legume crops. Therefore, studies on the root nodule development and senescence are of great significance in the efficient use of nitrogen sources and legume crop yield.

Legume nodule development is initiated by an exchange of chemical signals between legumes and soil rhizobia, and included lots of protein-protein interactions between legumes and soil rhizobia proteins (Zhang et al., 2018). The leguminous molecular events involved in nodulation and nodule development are quite well understood (Li Q. G. et al., 2013; Mergaert et al., 2020; Roy et al., 2020). The characteristics and the corresponding detection methods, the key leguminous genes currently known to be involved in the regulation of nodule senescence and various abiotic factors affecting nodule senescence were summarized in our previous review (Zhou et al., 2021). Except for legumes, the soil rhizobia also play critical roles in symbiotic nitrogen fixation progress (Lindström and Mousavi, 2020). Because Nod factors (NFs), surface polysaccharides and secreted proteins are three main determinants of host specificity in most rhizobia, the rhizobial genes that affect the biological synthesis of these signaling molecules were explored in nodulation, including nod, nif, fix genes, surface polysaccharides biosynthesis genes and secretion system-related genes (Wang et al., 2014; Li et al., 2020). Besides, it has been demonstrated that several rhizobial genes play key roles in nodule senescence. SmgshB is a critical gene for the synthesis of glutathione, which is a carbon source that can ensure the normal development of nodules and can suppress the early senescence of nodules (Cheng et al., 2017; Yang et al., 2020). lrpL-acdS genes, which encode an ACC deaminase enzyme, are associated with nitrogen fixation and delayed nodule senescence (Tittabutr et al., 2015). The sigma factors RpoH1 and RpoH2 are involved in different stress responses and negatively regulate nodule senescence (Martínez-Salazar et al., 2009). Anabaena variabilis flavodoxin, a protein involved in the response to oxidative stress, can significantly delay nodule senescence in alfalfa (Redondo et al., 2009). Therefore, both leguminous genes and rhizobial genes are critical to nodule development and senescence.

Soybean (Glycine max) is an important leguminous crop and has been grown worldwide for edible oil, food and feed material (Palander et al., 2005). Large amounts of N from the atmosphere fixed by soybean symbiotic nitrogen fixation can be widely available for growth of soybean and other intercropping and rotation crops (Biswas and Gresshoff, 2014). Branching stage, flowering stage, fruiting stage, pod stage and harvest stage are five important developmental stages for soybean cultivation and fertilizer application studies, and the nitrogen fixation characteristics of these stages have been well studied (Fehr et al., 1971). In our previous study, the expression of nodule genes at above-mentioned five stages was assessed quantitatively using RNA-Seq, we only analyzed the important different expression genes (DEGs) of soybean in nodule development and senescence, but not of Bradyrhizobium diazoefficiens 113-2 (Yuan et al., 2017).

B. diazoefficiens 113-2 was collected from paddy fields in Hengyang area of Hunan Province, China in 1972 by Xuejiang Zhang, and its genome shared a large synteny blocks and a high ANI value with that of B. diazoefficiens USDA110 (Li et al., 2020). In this report, we investigated the expression of B. diazoefficiens 113-2 genes in the nodule samples from above-mentioned five developmental stages of soybean and identified 164 important DEGs associated with nodule development and senescence. These DEGs included those encoding nod, nif, fix proteins and T3SS secretion system-related proteins, as well as proteins involved in nitrogen metabolism, ABC transporters and two-component system pathways. Our results firstly connected the action mode of rhizobial genes during nodule development with the developmental stages of soybean, which should provide new clues for the understanding of the genetic determinants of soil rhizobia in nodule development and senescence, and shed new light on the molecular mechanisms of nitrogen fixation at above-mentioned five developmental stages of soybean.

Results

Expression Analysis of the Genes of B. diazoefficiens 113-2

We previously performed RNA-Seq analysis for five different nodule samples from five important development stages of soybean (Yuan et al., 2017) and sequenced the genome of B. diazoefficiens 113-2 (Li et al., 2020). In this report, we want to investigate the expression of the genes of B. diazoefficiens 113-2 in nodule development and senescence by using the above-mentioned RNA-Seq data. The genome and gene mapping details of B. diazoefficiens 113-2 were shown in Supplementary Table 1, and the results of sequencing saturation, reads coverage and reads random analyses were displayed in Supplementary Figures 1–3. As shown in Table 1, 2234 (about 25.38%) genes of B. diazoefficiens 113-2 were not identified in these five nodule samples, and more than 30% of the genes (except for Pod stage_N) contained meaningless fpkm values (0∼1). For the rest genes, most of the fpkm values were between 1 and 1,000, and a very small number of the fpkm values were above 10,000 (Table 1). Besides, the numbers of genes with meaningful fpkm values (>1) at one or more nodule samples were shown in Figure 1A. There are 966 genes consistently expressed at five nodule samples, 1085 genes were consistently found at four nodule samples, and 1252, 1552, and 1771 genes expressed at three, two and only one nodule sample, respectively. The detail fpkm information was shown in Supplementary Table 2.

TABLE 1

fpkm value
>10,0001,000∼10,0001∼1,0000∼1NA
Branching stage_N16 (0.18%)223 (2.53%)2087 (23.71%)4241 (48.19%)2234 (25.38%)
Flowering stage_N18 (0.20%)181 (2.06%)2362 (26.84%)4006 (45.52%)2234 (25.38%)
Fruiting stage_N18 (0.20%)208 (2.36%)3293 (37.42%)3048 (34.63%)2234 (25.38%)
Pod stage_N18 (0.20%)210 (2.39%)5260 (59.77%)1079 (12.26%)2234 (25.38%)
Harvest stage_N25 (0.28%)190 (2.16%)3632 (41.27%)2720 (30.91%)2234 (25.38%)

The numbers of B. diazoefficiens 113-2 genes in different fpkm value groups in five nodule samples.

FIGURE 1

To study the correlation of the expressions of B. diazoefficiens 113-2 genes among the five nodule samples, we calculated the Pearson correlation coefficients of each two nodule samples based on the total gene expression levels (Figure 1B). The Pearson correlation coefficients between Branching stage_N and the other four nodule samples were 0.54∼0.63, indicating that the expression of the genes of B. diazoefficiens 113-2 in Branching stage_N was very different from other four nodule samples. The Pearson correlation coefficients between harvest stage_N and flowering stage_N, fruiting stage_N or pod stage_N were less than 0.9, suggesting that relative lower correlation between these nodule samples. While the Pearson correlation coefficients among flowering stage_N, fruiting stage_N and pod stage_N were more than 0.9, indicating that relative higher correlation were existed in these three nodule samples.

Different Expression Gene Identification Between Different Nodule Samples

To screen the important genes of B. diazoefficiens 113-2 involved in nodule development and senescence, a false discovery rate (FDR) ≤ 0.001 and | log2 ratio| ≥ 1 were used to identify the DEGs in the ten Groups. The numbers of up-regulated and down-regulated DEGs in these ten Groups were shown in Figure 2A, and the ID information on these DEGs was provided in Supplementary Table 3. The Group (Harvest vs. Branching) possessed the highest number of DEGs (36 up-regulated, 47 down-regulated), followed by Pod vs. Branching Group (33 up-regulated, 47 down-regulated) and Fruiting vs. Branching Group (34 up-regulated, 31 down-regulated). Besides, Harvest vs. Flowering Group (27 up-regulated, 25 down-regulated), Harvest vs. Pod Group (28 up-regulated, 13 down-regulated) and Harvest vs. Fruiting Group (19 up-regulated, 18 down-regulated) also had relative high number of DEGs. There were relatively low DEG numbers among the three Groups (Fruiting vs. Flowering, Pod vs. Flowering and Pod vs. Fruiting).

FIGURE 2

According to different reference developmental stages N, above-mentioned ten Groups were divided into four Classes (Figure 2A). The number of DEGs consistently found in all four Classes was 15 and the total number of DEGs identified during the soybean nodule development was 164 (Figure 2B).

DEG Annotation and Expression Pattern Analysis

The above-mentioned 164 DEGs encoded peptides with 41∼3105 amino acid residues, and their detail sequence information was shown in Supplementary Table 4. These DEGs were annotated to evaluate their potential functions. Among them, five DEGs-encoding nodulation proteins (NodG, NolG, NoeK, and NoeJ), nine DEGs-encoding Nif proteins (NifB, NifD, NifE, NifH, NifK, NifN, NifQ, NifW, and NifX), three DEGs-encoding Fix proteins (FixA, FixC, and FixK), three DEGs-encoding transcriptional regulator proteins and nine DEGs-encoding ABC transporter-related proteins were identified. Besides, 27 DEGs were involved in type-III secretion system (T3SS), eight DEGs located in membrane and many DEGs encoded various proteases (Supplementary Table 5).

To investigate the expression profile of the 164 DEGs during nodule development and senescence, we analyzed the RNA-Seq results of these DEGs at above-mentioned five nodule samples, and the result showed that these DEGs had various expression patterns (Figure 3). 64 DEGs reached their peaks at branching stage_N or flowering stage_N, and most of them were down-regulated during nodule development and senescence (Figure 3A), indicating that these DEGs mainly play roles in the early stages of nodule development. 45 DEGs reached their peaks at flowering stage_N or fruiting stage_N, and had relative low expression at branching stage_N or harvest stage_N (Figure 3B), suggesting that these DEGs mainly participate in nitrogen fixation progress. About 55 DEGs reached their peaks at harvest stage_N or pod stage_N, and most of them were up-regulated during nodule development and senescence (Figure 3C), indicating that these DEGs mainly play roles in nodule senescence.

FIGURE 3

DEG Functional Ontology and Kyoto Encyclopedia of Genes and Genomes Pathway Enrichment Analysis

To evaluate the potential functions of the 164 DEGs, we used Gene Ontology (an internationally standardized gene function classification system) to classify these DEGs. Among these DEGs, 62 genes had no functional GO term, and a total of 150 GO terms were assigned for the rest genes, while about 2/3 (101 out of 150) functional GO terms were specific (only for one gene) (Supplementary Table 6). The rest49 Go terms were listed in Figure 4 and were divided into three categories: cellular components, biological process and molecular function. The cellular components mainly included integral component of membrane and membrane. The main biological processes of the DEGs were oxidation-reduction process, nitrogen fixation, metabolic process, transport and regulation of transcription. The molecular functions associated with the DEGs mainly focused on oxidoreductase activity, catalytic activity, ATP binding, nucleic acid binding and DNA binding.

FIGURE 4

KEGG is the major public database for pathway enrichment analysis (Kanehisa et al., 2016). In this report, we used KEGG database to classify these DEGs. Among these DEGs, 74 genes had no classified KEGG pathways, and a total of 149 KEGG pathway subgroups were assigned for the rest genes, while about 76% (113 out of 149) pathway subgroups were specific (only for one gene) (Supplementary Table 7). The other 36 pathway subgroups associated with the DEGs were shown in Figure 5. Most of these pathways were metabolic-related pathways. The pathways with the greatest numbers were metabolic pathway, followed by microbial metabolism in diverse environments, biosynthesis of secondary metabolites and biosynthesis of antibiotics.

FIGURE 5

DEGs Associated With Nitrogen Metabolism, ABC Transporters and Two-Component System Pathways

To investigate the regulation of the DEGs of 113-2 during nodule development and senescence, we analyzed nitrogen metabolism (K00910), ABC transporters (K02010) and Two-Component system (K02020) pathways (obtained by KEGG)1 in more detail (Figure 6). Two KEGG gene sets in K00910, seven in K02010 and six in K02020 were differentially expressed during nodule development and senescence. One gene that matched K00910 pathway gene NirK (113-2GL002109) was up-regulated in Harvest vs. Branching, Harvest vs. Flowering, Harvest vs. Fruiting and Harvest vs. Pod, suggesting that this gene has relative higher expression in harvest stage nodules and plays key roles in nodule senescence. Among the three DEGs that matched k02010 pathway gene nifDKH, 113-2GL007881 has different expression patterns with 113-2GL007904 or 113-2GL007905 in Flowering vs. Branching, Fruiting vs. Branching, Pod vs. Branching and Harvest vs. Branching groups, indicating that nifDKH gene may play diverse roles in nodulation and nodule development. For K02010 pathway gene sets, three genes that matched modC, modA, gltI, and aatJ were down-regulated in all of the detected groups, the other gene that matched pstS was up-regulated in Pod vs. Branching and Pod vs. Flowering groups, and the rest two genes lptG and urtA were down-regulated then up-regulated during nodule development and senescence, suggesting that they may have roles in both nodulation and nodule senescence. Three K02020 pathway gene sets (pstS, gltI and aatJ) were the same as that of K02010 pathway, and among the rest three gene sets, one DEG that matched fixK (113-2GL006786) was up-regulated in all of the detected groups, one gene that matched ccoN (113-2GL006780) was up-regulated in most of the detected groups and only down-regulated in Fruiting vs. Flowering group, while the rest one DEG that matched atoB (113-2GL005802) was down-regulated in most of the detected groups and only up-regulated in Flowering vs. Branching group.

FIGURE 6

DEGs Encoding Nod and T3SS Secretion System-Related Proteins

In most rhizobia, NFs, surface polysaccharides and secreted proteins are the main determinants of host specificity (Putnoky et al., 1988; Lorkiewicz, 1997; Li et al., 2014). To examine whether these three signaling molecules also play roles during nodule development and senescence, we explored DEGs that required for the biological synthesis of these signaling molecules between different nodules from different developmental periods of soybean, including five Nod genes, eight nif genes, three fix genes and 27 T3SS-related genes (Figure 7). 16 DEGs (two Nod genes, seven nif genes, one fix gene and six T3SS-related genes) were down-regulated in all of the detected groups, suggesting that these genes mainly play roles in nodule development rather than nodule senescence. Ten DEGs (one Nod gene, one fix gene and eight T3SS-related genes) were up-regulated in all of the detected groups, and most of these genes (except for 113-2GL007701) have relative higher expression in pod stage nodules or harvest stage nodules, meaning that these genes may participate in regulating nodule senescence. Seven DEGs (113-2GL007881, 113-2GL007641, 113-2GL007564, 113-2GL007604, 113-2GL007623, 113-2GL007879, and 113-2GL007892) have highest expression in flowering stage nodules or fruiting stage nodules, reflecting that these genes may be critical to nitrogen fixation progress. Four DEGs (113-2GL001050, 113-2GL008064, 113-2GL002555, and 113-2GL004188) were down-regulated then up-regulated during nodule development and senescence, suggesting that these genes may play roles in both in the early stages of nodule development and nodule senescence. Five genes (113-2GL000406, 113-2GL005174, 113-2GL006780, 113-2GL007717, and 113-2GL008548) have relative high expression in flowering, fruiting or pod stage nodules, indicating that they mainly play roles in nitrogen fixation progress. The rest one gene (113-2GL007717) was only down-regulated in Harvest/Flowering group.

FIGURE 7

Orthologs of the DEGs

To further study the roles of these rhizobial DEGs during nodule development and senescence, we identified their orthologs in 113-2 itself, USDA110, Mesorhizobium loti MAFF303099, Mesorhizobium huakuii 7653R and Sinorhizobium meliloti 2011. Four homologous pairs (113-2GL000406/113-2GL001850, 113-2GL004736/113-2GL006371, 113-2GL004195/113-2GL002328, and 113-2GL003941/113-2GL001604) were existed in these 164 DEGs. 16 genes have no orthologs, and most of them have no COG functions (Supplementary Table 8), indicating that these DEGs were mainly strain-specific. 84 DEGs have only orthologs in USDA110, suggesting that these DEGs were mainly species-specific or host specific, and among these DEGs, about 13.1% (11 out of 84) have orthologs in 113-2 itself and more than one orthologs in USDA110 (Table 2). Among the rest 62 DEGs, 35 DEGs have one or more orthologs in all of the five rhizobia, 11 DEGs have orthologs in four rhizobia and 16 DEGs have orthologs in three rhizobial, reflecting that not all of these nodule development and senescence-related DEGs are core genes among these four rhizobial (Table 3). Among the core 35 DEGs, about 60% (21 out 35) were mainly expressed in the early stages of nodule development, six DEGs showed highest expression in the progress of nitrogen fixation, and eight DEGs may play critical roles in nodule senescence, indicating that these core genes have diverse roles during nodule development and senescence (Table 3). Among these four rhizobia, 113-2 or USDA110 and M. loti MAFF303099 can form determinate nodules with soybean and Lotus, respectively, and nine DEGs of 113-2 have orthologs only in USDA110 and M. loti MAFF303099 (Table 3). Besides, we selected four gene groups with more orthologs to perform phylogeny analyses, and the results revealed closer phylogenetic relationships between the two strains in the same genus (Supplementary Figures 4–7).

TABLE 2

DEGs113-2USDA110COG annotion
113-2GL000027bll0805CubicO group peptidase, beta-lactamase class C family
113-2GL000144blr0694Serine protease, subtilase family
113-2GL000369blr0497–
113-2GL000546bll0333Glucose dehydrogenase
113-2GL000547bll0332Cytochrome c, mono- and diheme variants
113-2GL000923113-2GL007956blr1693,bll8244–
113-2GL000935bll8232L-lysine 2,3-aminomutase (EF-P beta-lysylation pathway)
113-2GL000936bll8231Gamma-glutamyltranspeptidase
113-2GL000937bll8230–
113-2GL000966bll8195Retron-type reverse transcriptase
113-2GL001135bll7952DNA-binding beta-propeller fold protein YncE
113-2GL001326bll7787–
113-2GL001540113-2GL005055,113-2GL005658bll4384,bll7600,blr3849ABC-type amino acid transport/signal transduction system, periplasmic component/domain
113-2GL001603bll7540Fe-S oxidoreductase
113-2GL001604113-2GL003941bll7539,blr5311–
113-2GL001653bll7492Predicted DNA-binding protein, contains Ribbon-helix-helix (RHH) domain
113-2GL001765bll7395Uncharacterized protein
113-2GL001774bll7385Cellulose biosynthesis protein BcsQ
113-2GL001885blr7289Flavin-dependent oxidoreductase, luciferase family
113-2GL002261blr6953ABC-type molybdate transport system, ATPase component
113-2GL002263blr6951ABC-type molybdate transport system, periplasmic component
113-2GL002327blr6889DNA-binding transcriptional regulator, MarR family
113-2GL002328113-2GL004195, 113-2GL004285bll5076,bll4983,bll6888–
113-2GL002767bll6479ABC-type branched-chain amino acid transport system, periplasmic component
113-2GL002823bll6433–
113-2GL003062bll6129Uncharacterized protein
113-2GL004088bll5176–
113-2GL004268bsl5002–
113-2GL004412113-2GL004612, 113-2GL002719,blr4994,blr5174,bll6524,Opacity protein and related surface antigens
113-2GL004272, 113-2GL004615,blr4699,bll4867,blr4701,
113-2GL001425, 113-2GL004090blr7695
113-2GL004429113-2GL005811bll4853,blr3712Outer membrane protein assembly factor BamA
113-2GL004621bsr4694–
113-2GL004679blr4646CBS domain
113-2GL004737bll4594Fatty-acid desaturase
113-2GL005073bll4367Glyoxylase or a related metal-dependent hydrolase, beta-lactamase superfamily II
113-2GL005144bll4305Ribonuclease G or E
113-2GL005174bll4278–
113-2GL005281blr4186Glyoxylase or a related metal-dependent hydrolase, beta-lactamase superfamily II
113-2GL005303bsl4167–
113-2GL005372bll4106Lipopolysaccharide export LptBFGC system, permease protein LptF
113-2GL005840blr3683Chaperonin GroEL (HSP60 family)
113-2GL005842blr3681Predicted metal-dependent hydrolase, TIM-barrel fold
113-2GL005843blr3680–
113-2GL005844blr3679–
113-2GL005845blr3678Ferredoxin-NADP reductase
113-2GL005846blr3677–
113-2GL005848blr3675D-arabinose 1-dehydrogenase, Zn-dependent alcohol dehydrogenase family
113-2GL006036blr3477Aspartate/methionine/tyrosine aminotransferase
113-2GL006040blr3474Arginine/lysine/ornithine decarboxylase
113-2GL006052blr3464Uncharacterized protein
113-2GL006113blr3404Heme-degrading monooxygenase HmoA and related ABM domain proteins
113-2GL006583bll2956Sugar lactone lactonase YvrE
113-2GL007613113-2GL000934, 113-2GL007832,blr2165,blr7550,blr1829,Transposase
113-2GL007909, 113-2GL000934,blr1911,blr4868,blr1716,
113-2GL007752, 113-2GL007853,blr1702,bll4642,blr1740,
113-2GL007932, 113-2GL004410blr8233,bll8293,blr2076,
blr1807
113-2GL006652blr2887Uncharacterized protein
113-2GL006786113-2GL007581, 113-2GL006786,bll3466,bll2109,bll2757cAMP-binding domain of CRP or a regulatory subunit of cAMP-dependent protein kinases
113-2GL006050
113-2GL006982bll2590Nucleotide-binding universal stress protein, UspA family
113-2GL007163bll2460Outer membrane receptor for ferrienterochelin and colicins
113-2GL007170blr2455Isocitrate lyase
113-2GL007564blr2132–
113-2GL007565blr2131Lysine/ornithine N-monooxygenase
113-2GL007582blr2108Non-ribosomal peptide synthetase component F
113-2GL007583blr2106Mannose-6-phosphate isomerase, cupin superfamily
113-2GL007612blr2077O-acetylhomoserine/O-acetylserine sulfhydrylase, pyridoxal phosphate-dependent
113-2GL007619bll2067–
113-2GL007623bll2063Acyl-CoA dehydrogenase related to the alkylation response protein AidB
113-2GL007660bll2012Opacity protein and related surface antigens
113-2GL007665bll2003Carbohydrate-selective porin OprB
113-2GL007693blr1971Dipeptidyl aminopeptidase/acylaminoacyl peptidase
113-2GL007717bll1944Opacity protein and related surface antigens
113-2GL007787113-2GL000704bll0184,blr1879DnaJ-class molecular chaperone with C-terminal Zn finger domain
113-2GL007806113-2GL001935, 113-2GL001936blr7243,blr1853,blr7242Cytochrome P450
113-2GL007808blr1851Uncharacterized protein YjbI, contains pentapeptide repeats
113-2GL007874bll1777Alkyl hydroperoxide reductase subunit AhpC (peroxiredoxin)
113-2GL007880blr1770–
113-2GL007884113-2GL008289blr1311,bll1766Outer membrane protein W
113-2GL007893bll1754–
113-2GL007894blr1753Ferredoxin subunit of nitrite reductase or a ring-hydroxylating dioxygenase
113-2GL007900blr1748–
113-2GL007901blr1747Predicted Fe-Mo cluster-binding protein, NifX family
113-2GL007903blr1745Nitrogenase molybdenum-iron protein, alpha and beta chains
113-2GL008083blr1496Uncharacterized conserved protein YjbJ, UPF0337 family
113-2GL008134bll1447Superfamily II DNA and RNA helicase
113-2GL008317bll1285–
113-2GL008548blr10811,2-phenylacetyl-CoA epoxidase, catalytic subunit
113-2GL008607113-2GL006492bll1028, blr3042DNA-directed RNA polymerase specialized sigma subunit, sigma24 family

Analysis of species-specific or host specific DEGs.

TABLE 3

DEGsOrthologs of the DEGs
Expression in soybean symbiosis
113-2USDA110MAFF3030997653RSM2011(according to Figure 3)
113-2GL000054bll0779MAFF_RS22785MCHK_RS30795SM2011_RS15560Early stages of nodule development
113-2GL000055blr0778MAFF_RS28315MCHK_RS02430SM2011_RS06535Early stages of nodule development
113-2GL000126bll0710MAFF_RS20750MCHK_RS28780SM2011_RS15760Early stages of nodule development
113-2GL000160blr0678MAFF_RS19575MCHK_RS27375SM2011_RS15250Nodule senescence
113-2GL000269blr0573MAFF_RS17030MCHK_RS25150SM2011_RS31340Nodule senescence
113-2GL000406113-2GL001850bll7322,blr0462MAFF_RS17955MCHK_RS26065SM2011_RS30815Nodule senescence
113-2GL000429bll0442MAFF_RS16950MCHK_RS25060SM2011_RS29850Early stages of nodule development
113-2GL000514113-2GL002082bsr7117, blr0365MAFF_RS13260MCHK_RS21570SM2011_RS24865Early stages of nodule development
113-2GL000663113-2GL005801bll0225,blr3725MAFF_RS16145SM2011_RS31010Nodule senescence
113-2GL001050blr8272,bll8126SM2011_RS09965Early stages of nodule development
113-2GL001126113-2GL004687blr4637,blr7961MCHK_RS12390,SM2011_RS03045Nitrogen fixation
MCHK_RS18535
113-2GL001654blr7491MCHK_RS11600Nodule senescence
113-2GL001841blr7327MAFF_RS27500Nodule senescence
113-2GL001851blr7321MAFF_RS27480SM2011_RS04580Nodule senescence
113-2GL001863bll7310MAFF_RS27520MCHK_RS18460SM2011_RS04570,Nodule senescence
SM2011_RS01895
113-2GL002091bsl7109SM2011_RS33980Nitrogen fixation
113-2GL002109blr7089MCHK_RS17580SM2011_RS03405Nodule senescence
113-2GL002468bll6756MAFF_RS00585SM2011_RS21450Nodule senescence
113-2GL003340113-2GL001478bll7648,bll5866MAFF_RS18965MCHK_RS27010Nitrogen fixation
113-2GL003434blr5774MCHK_RS17320Nitrogen fixation
113-2GL003840bll5412MAFF_RS02150MCHK_RS09810SM2011_RS21190Early stages of nodule development
113-2GL003842bll5410MAFF_RS02160MCHK_RS09820SM2011_RS21200Early stages of nodule development
113-2GL003843bll5409MAFF_RS02165MCHK_RS09825SM2011_RS21205Early stages of nodule development
113-2GL003850bll5402MAFF_RS02200,SM2011_RS21230,Early stages of nodule development
MAFF_RS02105SM2011_RS21150
113-2GL003858bll5394MAFF_RS02240MCHK_RS09900SM2011_RS21270Early stages of nodule development
113-2GL003859bll5393MAFF_RS02245MCHK_RS09905SM2011_RS21275Early stages of nodule development
113-2GL003867bll5385MAFF_RS02285MCHK_RS09945SM2011_RS21315Early stages of nodule development
113-2GL004505bll4784MAFF_RS27120MCHK_RS01750SM2011_RS26345Nodule senescence
113-2GL004572113-2GL005733blr4723,blr3787MCHK_RS07395Nodule senescence
113-2GL004736113-2GL001679,bsl4595,bsr7468,MAFF_RS33280,MCHK_RS19875,SM2011_RS02010,Early stages of nodule development
113-2GL006371,bsl4595,bsl1445,MAFF_RS00310,MCHK_RS07180,SM2011_RS24070,
113-2GL008136bsr3154MAFF_RS11065,MCHK_RS19840,SM2011_RS00515,
MAFF_RS11510MCHK_RS32390SM2011_RS26245
MAFF_RS25580,MCHK_RS19440,SM2011_RS25125,
MAFF_RS11485MCHK_RS31800SM2011_RS24855
113-2GL005802bll0226,blr3724MAFF_RS16140MCHK_RS24370,SM2011_RS31015Nitrogen fixation
MCHK_RS02525
113-2GL005957blr3566MAFF_RS21725MCHK_RS29825Nitrogen fixation
113-2GL006512bll3022MAFF_RS15360MCHK_RS23560Nodule senescence
113-2GL006780blr2763MAFF_RS27090,MCHK_RS01720SM2011_RS02090,Nitrogen fixation
MAFF_RS26180SM2011_RS01660,
SM2011_RS03325
113-2GL007543blr2149MAFF_RS26000Early stages of nodule development
113-2GL007545blr2147MAFF_RS25985Early stages of nodule development
113-2GL007548blr2145MAFF_RS25975Nitrogen fixation
113-2GL007604113-2GL007555,bll3527,bll1900,MAFF_RS25425MCHK_RS31900SM2011_RS30325Nitrogen fixation
113-2GL003457,bll5000,blr1706
113-2GL008356,bll1997
113-2GL007792,
113-2GL004069,
113-2GL000973,
113-2GL007670,
113-2GL007896,
113-2GL007767,
113-2GL007630,
113-2GL007942,
113-2GL004713,
113-2GL007823
113-2GL007618blr2068MAFF_RS25050Nodule senescence
113-2GL007625113-2GL002237,blr5226,blr6978,MAFF_RS33345,MCHK_RS17610,SM2011_RS18355,Early stages of nodule development
113-2GL007625,blr5625,blr5226,MAFF_RS09905,MCHK_RS07245SM2011_RS00340,
113-2GL001613,bll2060,bsr7532MAFF_RS10535,SM2011_RS02030,
113-2GL004032,MAFF_RS36295,SM2011_RS11615
113-2GL003601MAFF_RS23815
113-2GL007626113-2GL001611,blr5626,bll2059,MAFF_RS23810,MCHK_RS07240SM2011_RS11610,Early stages of nodule development
113-2GL004031,blr6979,blr5227,MAFF_RS36300,SM2011_RS20400,
113-2GL003600,blr4635,blr7533MAFF_RS33340,SM2011_RS18350,
113-2GL002236,MAFF_RS10540SM2011_RS02025
113-2GL004689
113-2GL007641blr2038MAFF_RS24010MCHK_RS32450SM2011_RS02275Nitrogen fixation
113-2GL007643113-2GL006940bll2623,blr2036MAFF_RS13735MCHK_RS22150SM2011_RS04845Nitrogen fixation
113-2GL007877blr1774MAFF_RS24000MCHK_RS32460SM2011_RS02265Early stages of nodule development
113-2GL007879blr1771MAFF_RS24015Nitrogen fixation
113-2GL007881blr1769MAFF_RS24195MCHK_RS32325SM2011_RS02280Nitrogen fixation
113-2GL007888blr1759MAFF_RS23985MCHK_RS32485SM2011_RS02250Early stages of nodule development
113-2GL007892blr1755MAFF_RS24355Nitrogen fixation
113-2GL007898bsr1750MAFF_RS24040,Early stages of nodule development
MAFF_RS24300
113-2GL007899bsr1749MAFF_RS24035,Early stages of nodule development
MAFF_RS24295
113-2GL007902blr1746MAFF_RS24215MCHK_RS32305SM2011_RS02420Early stages of nodule development
113-2GL007904blr1744MAFF_RS24205MCHK_RS32315SM2011_RS02290Early stages of nodule development
113-2GL007905blr1743MAFF_RS24200MCHK_RS32320SM2011_RS02285Early stages of nodule development
113-2GL007964113-2GL005342blr1686,blr4134SM2011_RS08430Nitrogen fixation
113-2GL008056bll1523MAFF_RS15765MCHK_RS33030,SM2011_RS28455Nodule senescence
MCHK_RS23965
113-2GL008064blr1516SM2011_RS15300Early stages of nodule development
113-2GL008106bsl1473MAFF_RS37185SM2011_RS11135Nitrogen fixation
113-2GL008132113-2GL003149blr1448,blr6053MAFF_RS04190MCHK_RS11590SM2011_RS14035Nodule senescence
113-2GL008135bsl1446MAFF_RS26365MCHK_RS00670SM2011_RS17430Early stages of nodule development
113-2GL008535blr1091MAFF_RS15660MCHK_RS23865Nodule senescence

Orthologs of the selected DEGs in 113-2, USDA110, M. loti MAFF303099, M. huakuii 7653R, and S. meliloti 2011.

Our previous study has indicated that 113-2 has an extremely close phylogenetic relationship to USDA110 (Li et al., 2020). An interspecies protein interactome was constructed and 290 USDA110 proteins associated with soybean symbiosis were provided (Zhang et al., 2018). In this report, about 203 USDA110 proteins were searched from the orthologs of the DEGs (Tables 2, 3). We compared these two USDA110 protein groups and discovered six same proteins (bll7322, blr1971, blr1091, blr1516, blr0462, and bll7600), which were detected to be expressed in rhizoidal of root nodules during symbiosis in previous studies (Pessi et al., 2007; Delmotte et al., 2010; ÄŒuklina et al., 2016; Zhang et al., 2018). Then we predict the interaction proteins of the six DEGs (orthologs of the six USDA110 proteins) according to the resources from the constructed protein interactome (Zhang et al., 2018), and the results were showed in Table 4. These rhizoidal proteins have more than two predicted interacting proteins in soybean, especially for bll7322 or 113-2GL000406 or 113-2GL001850, and blr1971 or 113-2GL007693; they have 21 and 25 predicted interacting proteins in soybean, respectively.

TABLE 4

DEGsOrthologs in USDA110References for orthologsPotencial interaction proteins
SoybeanUSDA110113-2
113-2GL000406bll7322Zhang et al., 2018Glyma.01G003800bll7785113-2GL002891
113-2GL001850Glyma.02G017800113-2GL001328
Glyma.02G147200
Glyma.02G222300
Glyma.03G164400
Glyma.03G251900
Glyma.04G080800
Glyma.04G082000
Glyma.04G194700
Glyma.04G226700
Glyma.06G081900
Glyma.06G082500
Glyma.06G171200
Glyma.09G273200
Glyma.10G018200
Glyma.12G045400
Glyma.12G167500
Glyma.16G043900
Glyma.18G284300
Glyma.19G222600
Glyma.20G249300
blr0462Pessi et al., 2007Glyma.08G152600
Delmotte et al., 2010Glyma.15G272000
ÄŒuklina et al., 2016
Zhang et al., 2018
113-2GL007693blr1971Pessi et al., 2007Glyma.01G118900blr1971113-2GL007693
Delmotte et al., 2010Glyma.02G155300bll5750113-2GL003467
ÄŒuklina et al., 2016Glyma.02G208700
Zhang et al., 2018Glyma.03G004900
Glyma.03G195400
Glyma.04G007400
Glyma.04G081100
Glyma.05G204000
Glyma.07G079700
Glyma.07G173100
Glyma.08G067400
Glyma.09G278600
Glyma.10G280000
Glyma.11G102100
Glyma.11G237400
Glyma.14G176900
Glyma.14G206300
Glyma.15G058400
Glyma.16G117600
Glyma.18G056100
Glyma.18G060600
Glyma.18G210100
Glyma.19G173400
Glyma.19G205200
Glyma.20G005900
113-2GL008535blr1091Pessi et al., 2007Glyma.02G294900blr1148113-2GL008473
Delmotte et al., 2010Glyma.14G018700blr3906113-2GL005583
ÄŒuklina et al., 2016Glyma.14G161800blr4307113-2GL005142
Zhang et al., 2018Glyma.14G162000bll4945113-2GL004326
113-2GL008064blr1516Pessi et al., 2007Glyma.05G033800blr1516113-2GL008064
Delmotte et al., 2010Glyma.14G145000bll3640113-2GL005884
ÄŒuklina et al., 2016blr4473113-2GL004883
Zhang et al., 2018bll4736113-2GL004557
bll5710113-2GL003508
113-2GL001540bll7600Delmotte et al., 2010Glyma.17G228800bll7600113-2GL001540
ÄŒuklina et al., 2016Glyma.18G277300
Zhang et al., 2018

The prediction of the interaction proteins of the selected DEGs.

Discussion

Symbiotic nitrogen fixation system is of great significance in agriculture and ecology, and nodule development and senescence directly affects nitrogen fixation efficiency. Nodule development and senescence involves a programmed series of molecular events that are carried out synchronously by both legumes and rhizobia (Vasse et al., 1990; Cermola et al., 2000; Oldroyd et al., 2011; Roux et al., 2014; Alemneh et al., 2020). The nitrogen fixation characteristics of five soybean developmental stages (Branching stage, flowering stage, fruiting stage, pod stage and harvest stage) have been well studied (Fehr et al., 1971), the expression of nodule genes at these developmental stages was assessed quantitatively using RNA-Seq, and the important DEGs of soybean that associated with nodule development and senescence were identified in our previous study (Yuan et al., 2017). In the present study, we firstly investigated the expression of rhizobial genes in the nodule samples from above-mentioned five developmental stages of soybean by using our previous RNA-Seq data (Yuan et al., 2017)and identified 164 important rhizobial DEGs associated with nodule development and senescence.

The Gene Expression of B. diazoefficiens 113-2 at Different Nodule Samples From Different Soybean Developmental Stages

Previous comparative transcriptomic analyses have explored the gene expression patterns of rhizobia during chemoautotrophic growth (Franck et al., 2008) or symbiotic non-growth (Vercruysse et al., 2011), under different stresses (Ampe et al., 2003; Chang et al., 2007; Cytryn et al., 2007), in determinate nodules and indeterminate nodules (Li Y. et al., 2013), and in different zones of the indeterminate nodule (Roux et al., 2014). For the bacteroids in soybean symbiosis, several transcriptomic analyses of bacteroids focused on symbiotic gene region (Hauser et al., 2006) and single time point during optimal symbiotic functioning (Chang et al., 2007; Pessi et al., 2007; Franck et al., 2014). It is worth noting that a multiple time-point microarray study in B. diazoefficiens bacteroids isolated from soybean nodules from nodulation to nodule senescence, which revealed a shift gene expression patterns of B. diazoefficiens during symbiotic process, especially in nodule senescence (Franck et al., 2018). In this report, we also performed a multiple time-point transcriptomic ananlysis, the difference is that the RNA samples were collected from different nodule samples from five important developmental stages of soybean (Yuan et al., 2017). We firstly investigated the expression of rhizobial genes in symbiosis at above-mentioned five developmental stages of soybean, which should shed new light on the molecular mechanisms of nitrogen fixation during the growth and development of soybean.

In the five soybean nodule samples, about 25.38% (2234 out of 8801) of B. diazoefficiens 113-2 genes were not detected, indicating that these genes may not participate in symbiotic progress. The pod stage nodules possessed the highest number of genes with meaningful fpkm values (>1), suggesting the initiation of a series of new processes, which might be associated with nodule senescence. The correlation analysis of the gene expressions in the five nodule samples showed that the expression of B. diazoefficiens 113-2 genes at Branching stage or harvest stage has relative lower correlation with the other soybean development stages, but high correlations were existed among flowering stage, fruiting stage and pod stage (Figure 1B), indicating that similar gene expression patterns of B. diazoefficiens 113-2 were found during optimal symbiotic functioning, but different expression patterns were existed among early nodule development, nitrogen fixation progress and nodule senescence. These results also revealed a shift gene expression pattern of B. diazoefficiens 113-2 in symbiosis during the growth and development of soybean.

DEGs of B. diazoefficiens 113-2 Associated With Nodule Development and Senescence

Previous study showed that > 800 USDA110 genes were up- or down-regulated at each time point of soybean nodule development (Franck et al., 2018), while in this report, the highest number of DEGs at the ten groups was 83 (Figure 2A). The major reason may be that we used our previous soybean nodule RNA-Seq data, in which the RNA samples were mixed with both soybean RNA and B. diazoefficiens 113-2 RNA (Yuan et al., 2017). Our main purpose is to screen the important DEGs that can interact with host and are mainly involved in nodule development and senescence. The gene expression of soil rhizobia inside of the host was different from free live state (Vercruysse et al., 2011), so we did not isolate B. diazoefficiens 113-2 from the nodule samples.

Lots of important rhizobial genes have been showed to participate in nodule development and senescence (Voroshilova et al., 2001; Tittabutr et al., 2015; Franck et al., 2018). Recently, we sequenced the genome of B. diazoefficiens 113-2 (Li et al., 2020) and identified 164 DEGs of B. diazoefficiens 113-2 during the soybean nodule development in this study. Among these DEGs, 49 DEGs had both no functional GO term and no classified KEGG pathway (Supplementary Tables 6, 7), and about 13 DEGs were not clearly annotated (Supplementary Table 5), meaning that we identified many novel and uncharacterized genes associated with nodule development and senescence. About 62.2% (102 out of 164) of the DEGs have no orthologs in M. loti MAFF303099, M. huakuii 7653R and S. meliloti2011, and only 35 DEGs have orthologs in all of the three rhizobia (Table 3), indicating that most of the identified genes associated with nodule development and senescence are not core genes among the rhizobia genomes. Besides, B. diazoefficiens 113-2 and M. loti MAFF303099 can form determinate nodules with soybean and Lotus, respectively, but only nine DEGs of 113-2 have orthologs only in M. loti MAFF303099 (Table 3), suggesting that the nodule development and senescence in determinate nodules or indeterminate nodules was also mainly determined by host legume plants, which is similar to our previous study (Li et al., 2020).

The Potential Roles of Host Specificity-Related DEGs in Nodule Development and Senescence

In most rhizobia, the genes that affect the biosynthesis of nodulation factors, secretion system and surface polysaccharides are often play critical roles in host specificity (Horvath et al., 1986; Philip-Hollingsworth et al., 1989; Wang et al., 2014). In this report, five Nod genes, eight nif genes, three fix genes and 27 T3SS-related genes (no genes associated with surface polysaccharides biosynthesis) were identified (Figure 7). The Nod genes primarily played roles in nodulation (Moulin et al., 2004), and among the identified five Nod genes, only Noek (113-2GL001050) and NolG (113-2GL008064) may play roles in the early state of nodule development. NoeJ (113-2GL007583) and NolG (113-2GL007643) may participate in nitrogen fixation progress and NodG (113-2000663) may have roles in nodule senescence. The nif and fix genes are crucial for nitrogen fixation progress (Fischer and Hennecke, 1984; Hennecke, 1990). NifH, NifD, NifK, NifE, and NifN proteins are the core components of nitrogenase (Raymond et al., 2004), and the FixABCX genes were found to encode a membrane complex participating in electron transfer to nitrogenase (Earl et al., 1987; Edgren and Nordlund, 2004). The two main regulatory cascades that can be found in rhizobia are the RpoN−NifA and the oxygen−responsive two−component FixL−FixJ system, together with FixK (Lindström and Mousavi, 2019), which was found to regulate positively and negatively nitrogen fixation genes in Rhizobium meliloti (Batut et al., 1989). Among the identified eight Nif genes, only NifH (113-2GL007881) and NifW (113-2GL007879) reached the expression peaks at flowering stage_N or fruiting stage_N, the rest six nif genes NifB (113-2GL007888), NifD/K (113-2GL007904, 113-2GL007095), NifE (113-2-GL007903), NifN (113-2GL007902), NifQ (113-2007880), and NifX (113-2GL007901) may mainly involved in regulating the early nodule development. For the three fix genes, FixA (113-2GL007641) and FixC (113-2GL007877) had relative higher expression at branching stage_N and flower stage_N, implying that they may mainly play roles in nitrogen fixation progress, while FixK (113-2GL006786) may mainly have roles in nodule senescence. The rhizobia T3SS and its secreted effectors has been reported to modulate nodulation and host range (Okazaki et al., 2013). For the 27 T3SS-related genes, only seven DEGs had relative higher expression at branching stage_N, and the rest 20 genes may participate in nitrogen fixation progress or nodule senescence, meaning that the T3SS-related genes have roles in the whole nodule development and senescence rather than only in nodulation, which is similar to our previous study (Yuan et al., 2017). Together, these results revealed that the identified host specificity-related DEGs have diversified roles in symbiosis, meaning that host specific regulation maybe existed in the whole progress of nodule development and senescence, not only in nodulation.

In summary, we firstly combined the expression characteristics of rhizobial genes in symbiosis with the growth and development of soybean, and identified 164 important rhizobial DEGs associated with nodule development and senescence. Among these DEGs, many are firstly identified and showed to associate with nodule development and senescence. Besides, we focused on the DEGs encoding nod, nif, fix proteins and T3SS secretion system-related proteins, as well as proteins involved in nitrogen metabolism, ABC transporters and two-component system pathways. Our results supply valuable basises for studying the genetic determinants of soil rhizobia in symbiosis and the creation of high efficiency soybean cultivation technology.

Materials and Methods

Genome and Gene Mapping

The clean reads for the five nodule samples at five developmental stages of soybean (branching stage, flowering stage, fruiting stage, pod stage and harvest stage) were obtained from our previous study (Yuan et al., 2017). The raw sequence reads have been submitted to NCBI under the assigned accession number PRJNA765164, and the BioSample accessions included SAMN21545854, SAMN21545855, SAMN21545856, SAMN21545857, and SAMN21545858. The information of the reference genome (B. diazoefficiens 113-2) was previously described (Li et al., 2020). We map clean reads to the reference genome using HISAT (Cock et al., 2010), and map clean reads to reference transcripts using Bowtie2 (Kim et al., 2015). The mapping details are shown in Supplementary Table 1.

Gene Expression Analysis and Correlation Analysis Between Samples

We calculate the gene expression level for each nodule sample with RSEM (Li and Dewey, 2011). Besides, in order to reflect the gene expression correlation between the nodule samples, we calculated the Pearson correlation coefficients for all gene expression levels between each two samples using cor function in R software, and reflected these coefficients in the form of heat maps. Heatmap was performed using TBtools software (Chen et al., 2020).

Screening of DEGs Between Different Nodule Samples

To screen the DEGs between nodule samples, in combination with the false discovery rata (FDR) ≤ 0.001 and | log2 ratio| ≥ 1, the fold changes in gene expression were used to evaluate the significance of gene expression differences in the ten Groups which are classified according to our previous study (Yuan et al., 2017). Data of the DEGs in the ten Groups between five nodule samples are provided in Supplementary Table 3.

Gene Ontology Functional and Kyoto Encyclopedia of Genes and Genomes Pathway Analyses of DEGs

Gene Ontology (GO) is a database created by the Gene Ontology Consortium and is a classification system for genes’ biological functions. The method of GO functional enrichment analysis first maps all DEGs to terms in the GO database2 and the computing method as previously described (Yuan et al., 2017). Kyoto Encyclopedia of Genes and Genomes (KEGG) is a public pathway-related database,3 and the computing method also were described in our previous study (Yuan et al., 2017).

Orthologs Analysis of the DEGs

Firstly, we performed Core/Pan gene analysis for B. diazoefficiens 113-2, M. loti MAFF303099, M. huakuii 7653R, and S. meliloti 2011. We used the CD-HIT 4.664 rapid clustering of similar proteins software (Edgar, 2004)with a threshold of 50% pairwise identity and 0.7 length difference cut off in amino acid to cluster the Core/Pan genes of these four strains, and obtained the pan gene pool. Then we used the identified DEGs as queries to search the orthologs in the pan gene pool.

Phylogenetic Analysis

Phylogenetic trees were conducted in MEGA X (Kumar et al., 2018). The evolutionary history was inferred using the Neighbor-Joining method (Saitou and Nei, 1987). The percentage of replicate trees in which the associated taxa clustered together in the bootstrap test (1000 replicates) are shown next to the branches.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Statements

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

Author contributions

SY, HC, and XZ designed this work. SY wrote the manuscript. SY, SZ, and YF performed most of the analysis, figures, and tables. CZ, YH, ZS, SC, WG, HY, ZY, and DQ contributed substantially to the completion of this work. All authors contributed to the article and approved the submitted version.

Funding

This research was funded by the National Natural Science Foundation of China (Grant No. 32071964), the National Key Research and Development Program of China (Grant No. 2020YFD1000903-10), and the Basic Scientific Research Service Fee Special of the Central Scientific Research Institute (Grant No. 1610172018001).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2021.754837/full#supplementary-material

References

  • 1

    AlemnehA. A.ZhouY.RyderM. H.DentonM. D. (2020). Mechanisms in plant growth-promoting rhizobacteria that enhance legume-rhizobial symbioses.J. Appl. Microbiol.1291133–1156. 10.1111/jam.14754

  • 2

    AmpeF.KissE.SabourdyF.BatutJ. (2003). Transcriptome analysis of Sinorhizobium meliloti during symbiosis.Genome Biol.4:R15. 10.1186/gb-2003-4-2-r15

  • 3

    BatutJ.Daveran-MingotM. L.DavidM.JacobsJ.GarneroneA. M.KahnD. (1989). fixK, a gene homologous with fnr and crp from Escherichia coli, regulates nitrogen fixation genes both positively and negatively in Rhizobium meliloti.EMBO J.81276–1286. 10.1002/j.1460-2075.1989.tb03502.x

  • 4

    BiswasB.GresshoffP. M. (2014). The role of symbiotic nitrogen fixation in sustainable production of biofuels.Int. J. Mol. Sci.157380–7397. 10.3390/ijms15057380

  • 5

    CermolaM.FedorovaE.TatéR.RiccioA.FavreR.PatriarcaE. J. (2000). Nodule invasion and symbiosome differentiation during Rhizobium etli-Phaseolus vulgaris symbiosis.Mol. Plant Microbe Interact.13733–741. 10.1094/MPMI.2000.13.7.733

  • 6

    ChangW. S.FranckW. L.CytrynE.JeongS.JoshiT.EmerichD. W.et al (2007). An oligonucleotide microarray resource for transcriptional profiling of Bradyrhizobium japonicum.Mol. Plant Microbe Interact.201298–1307. 10.1094/MPMI-20-10-1298

  • 7

    ChenC.ChenH.ZhangY.ThomasH. R.FrankM. H.HeY.et al (2020). TBtools: an Integrative Toolkit Developed for Interactive Analyses of Big Biological Data.Mol. Plant131194–1202. 10.1016/j.molp.2020.06.009

  • 8

    ChengG.KarunakaranR.EastA. K.Munoz-AzcarateO.PooleP. S. (2017). Glutathione affects the transport activity of Rhizobium leguminosarum 3841 and is essential for efficient nodulation.FEMS Microbiol. Lett.364:fnx045. 10.1093/femsle/fnx045

  • 9

    CockP. J.FieldsC. J.GotoN.HeuerM. L.RiceP. M. (2010). The Sanger FASTQ file format for sequences with quality scores, and the Solexa/Illumina FASTQ variants.Nucleic Acids Res.381767–1771. 10.1093/nar/gkp1137

  • 10

    ČuklinaJ.HahnJ.ImakaevM.OmasitsU.FörstnerK. U.LjubimovN.et al (2016). Genome-wide transcription start site mapping of Bradyrhizobium japonicum grown free-living or in symbiosis - a rich resource to identify new transcripts, proteins and to study gene regulation.BMC Genomics17:302. 10.1186/s12864-016-2602-9

  • 11

    CytrynE. J.SangurdekarD. P.StreeterJ. G.FranckW. L.ChangW. S.StaceyG.et al (2007). Transcriptional and physiological responses of Bradyrhizobium japonicum to desiccation-induced stress.J. Bacteriol.1896751–6762. 10.1128/JB.00533-07

  • 12

    DelmotteN.AhrensC. H.KniefC.QeliE.KochM.FischerH. M.et al (2010). An integrated proteomics and transcriptomics reference data set provides new insights into the Bradyrhizobium japonicum bacteroid metabolism in soybean root nodules.Proteomics101391–1400. 10.1002/pmic.200900710

  • 13

    EarlC. D.RonsonC. W.AusubelF. M. (1987). Genetic and structural analysis of the Rhizobium meliloti fixA, fixB, fixC, and fixX genes.J. Bacteriol.1691127–1136. 10.1128/jb.169.3.1127-1136.1987

  • 14

    EdgarR. C. (2004). MUSCLE: a multiple sequence alignment method with reduced time and space complexity.BMC Bioinformatics5:113. 10.1186/1471-2105-5-113

  • 15

    EdgrenT.NordlundS. (2004). The fixABCX genes in Rhodospirillum rubrum encode a putative membrane complex participating in electron transfer to nitrogenase.J. Bacteriol.1862052–2060. 10.1128/JB.186.7.2052-2060.2004

  • 16

    FehrW. R.CavinessC. E.BurmoodD. T.PenningtonJ. S. (1971). Stage of Development Descriptions for Soybeans, Glycine Max (L.) Merrill1.Crop Sci.11929–931. 10.2135/cropsci1971.0011183X001100060051x

  • 17

    FergusonB. J.IndrasumunarA.HayashiS.LinM. H.LinY. H.ReidD. E.et al (2010). Molecular analysis of legume nodule development and autoregulation.J. Integr. Plant Biol.5261–76. 10.1111/j.1744-7909.2010.00899.x

  • 18

    FischerH. M.HenneckeH. (1984). Linkage map of the Rhizobium japonicum nifl-I and nitDK operons encoding the polypeptides of the nitrogenase enzyme complex.Mol. Gen. Genet.196537–540. 10.1007/BF00436206

  • 19

    FranckS.FranckW. L.BirkeS. R.ChangW.-S.SangurdekarD. P.CytrynE.et al (2014). Comparative transcriptomic analysis of symbiotic Bradyrhizobium japonicum.Symbiosis63123–135. 10.1007/s13199-014-0294-y

  • 20

    FranckS.StrodtmanK. N.QiuJ.EmerichD. W. (2018). Transcriptomic Characterization of Bradyrhizobium diazoefficiens Bacteroids Reveals a Post-Symbiotic, Hemibiotrophic-Like Lifestyle of the Bacteria within Senescing Soybean Nodules.Int. J. Mol. Sci.19:3918. 10.3390/ijms19123918

  • 21

    FranckW. L.ChangW. S.QiuJ.SugawaraM.SadowskyM. J.SmithS. A.et al (2008). Whole-genome transcriptional profiling of Bradyrhizobium japonicum during chemoautotrophic growth.J. Bacteriol.1906697–6705. 10.1128/JB.00543-08

  • 22

    HauserF.LindemannA.VuilleumierS.PatrignaniA.SchlapbachR.FischerH. M.et al (2006). Design and validation of a partial-genome microarray for transcriptional profiling of the Bradyrhizobium japonicum symbiotic gene region.Mol. Genet. Genomics27555–67. 10.1007/s00438-005-0059-7

  • 23

    HenneckeH. (1990). Nitrogen fixation genes involved in the Bradyrhizobium japonicum-soybean symbiosis.FEBS Lett.268422–426. 10.1016/0014-5793(90)81297-2

  • 24

    HirschA. M. (1992). Tansley Review No. 40. Developmental Biology of Legume Nodulation.New Phytol.122211–237. 10.1111/j.1469-8137.1992.tb04227.x

  • 25

    HorvathB.KondorosiE.JohnM.SchmidtJ.TörökI.GyörgypalZ.et al (1986). Organization, structure and symbiotic function of Rhizobium meliloti nodulation genes determining host specificity for alfalfa.Cell46335–343. 10.1016/0092-8674(86)90654-9

  • 26

    KanehisaM.SatoY.KawashimaM.FurumichiM.TanabeM. (2016). KEGG as a reference resource for gene and protein annotation.Nucleic Acids Res.44D457–D462. 10.1093/nar/gkv1070

  • 27

    KarrD. B.EmerichD. W. (1988). Uniformity of the microsymbiont population from soybean nodules with respect to buoyant density.Plant Physiol.86693–699. 10.1104/pp.86.3.693

  • 28

    KarrD. B.SuzukiF.WatersJ. K.EmerichD. W. (1990). Further Evidence for the Uniformity of the Microsymbiont Population from Soybean Nodules.J. Plant Physiol.136659–663. 10.1016/S0176-1617(11)81340-4

  • 29

    KimD.LangmeadB.SalzbergS. L. (2015). HISAT: a fast spliced aligner with low memory requirements.Nat. Methods12357–360. 10.1038/nmeth.3317

  • 30

    KumarS.StecherG.LiM.KnyazC.TamuraK. (2018). MEGA X: molecular Evolutionary Genetics Analysis across Computing Platforms.Mol. Biol. Evol.351547–1549. 10.1093/molbev/msy096

  • 31

    LiB.DeweyC. N. (2011). RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome.BMC Bioinformatics12:323. 10.1186/1471-2105-12-323

  • 32

    LiQ. G.ZhangL.LiC.DunwellJ. M.ZhangY. M. (2013). Comparative genomics suggests that an ancestral polyploidy event leads to enhanced root nodule symbiosis in the Papilionoideae.Mol. Biol. Evol.302602–2611. 10.1093/molbev/mst152

  • 33

    LiY.TianC. F.ChenW. F.WangL.SuiX. H.ChenW. X. (2013). High-resolution transcriptomic analyses of Sinorhizobium sp. NGR234 bacteroids in determinate nodules of Vigna unguiculata and indeterminate nodules of Leucaena leucocephala.PLoS One8:e70531. 10.1371/journal.pone.0070531

  • 34

    LiR.FengY.ChenH.ZhangC.HuangY.ChenL.et al (2020). Whole-Genome Sequencing of Bradyrhizobium diazoefficiens 113-2 and Comparative Genomic Analysis Provide Molecular Insights Into Species Specificity and Host Specificity.Front. Microbiol.11:576800. 10.3389/fmicb.2020.576800

  • 35

    LiX.DengZ.LiuZ.YanY.WangT.XieJ.et al (2014). The genome of Paenibacillus sabinae T27 provides insight into evolution, organization and functional elucidation of nif and nif-like genes.BMC Genomics15:723. 10.1186/1471-2164-15-723

  • 36

    LindströmK.MousaviS. A. (2019). Effectiveness of nitrogen fixation in rhizobia.Microb. Biotechnol.131314–1335. 10.1111/1751-7915.13517

  • 37

    LindströmK.MousaviS. A. (2020). Effectiveness of nitrogen fixation in rhizobia.Microb. Biotechnol.131314–1335.

  • 38

    LorkiewiczZ. (1997). Nodulation genes in the Rhizobium–plant signal exchange.Acta Biochim. Pol.441–12. 10.18388/abp.1997_4434

  • 39

    Martínez-SalazarJ. M.Sandoval-CalderónM.GuoX.Castillo-RamírezS.ReyesA.LozaM. G.et al (2009). The Rhizobium etli RpoH1 and RpoH2 sigma factors are involved in different stress responses.Microbiology155386–397. 10.1099/mic.0.021428-0

  • 40

    MergaertP.KeresztA.KondorosiE. (2020). Gene Expression in Nitrogen-Fixing Symbiotic Nodule Cells in Medicago truncatula and Other Nodulating Plants.Plant Cell3242–68. 10.1105/tpc.19.00494

  • 41

    MoulinL.BénaG.Boivin-MassonC.StêpkowskiT. (2004). Phylogenetic analyses of symbiotic nodulation genes support vertical and lateral gene co-transfer within the Bradyrhizobium genus.Mol. Phylogenet. Evol.30720–732. 10.1016/S1055-7903(03)00255-0

  • 42

    OkazakiS.KanekoT.SatoS.SaekiK. (2013). Hijacking of leguminous nodulation signaling by the rhizobial type III secretion system.Proc. Natl. Acad. Sci. U. S. A.11017131–17136. 10.1073/pnas.1302360110

  • 43

    OldroydG. E.MurrayJ. D.PooleP. S.DownieJ. A. (2011). The rules of engagement in the legume-rhizobial symbiosis.Annu. Rev. Genet.45119–144. 10.1146/annurev-genet-110410-132549

  • 44

    PalanderS.NäsiM.ValkonenE. (2005). Effect of age of growing turkeys on nutrient digestibility and energy value of cereal-soybean-based diets.Arch. Anim. Nutr.59139–147. 10.1080/17450390512331387945

  • 45

    PessiG.AhrensC. H.RehrauerH.LindemannA.HauserF.FischerH. M.et al (2007). Genome-wide transcript analysis of Bradyrhizobium japonicum bacteroids in soybean root nodules.Mol. Plant Microbe Interact.201353–1363. 10.1094/MPMI-20-11-1353

  • 46

    Philip-HollingsworthS.HollingsworthR. I.DazzoF. B.DjordjevicM. A.RolfeB. G. (1989). The effect of interspecies transfer of Rhizobium host-specific nodulation genes on acidic polysaccharide structure and in situ binding by host lectin.J. Biol. Chem.2645710–5714. 10.1016/S0021-9258(18)83607-9

  • 47

    PutnokyP.GrosskopfE.HaD. T.KissG. B.KondorosiA. (1988). Rhizobium fix genes mediate at least two communication steps in symbiotic nodule development.J. Cell Biol.106597–607. 10.1083/jcb.106.3.597

  • 48

    RaymondJ.SiefertJ. L.StaplesC. R.BlankenshipR. E. (2004). The Natural History of Nitrogen Fixation.Mol. Biol. Evol.21541–554. 10.1093/molbev/msh047

  • 49

    RedondoF. J.De La PeñaT. C.MorcilloC. N.LucasM. M.PueyoJ. J. (2009). Overexpression of flavodoxin in bacteroids induces changes in antioxidant metabolism leading to delayed senescence and starch accumulation in alfalfa root nodules.Plant Physiol.1491166–1178. 10.1104/pp.108.129601

  • 50

    RouxB.RoddeN.JardinaudM. F.TimmersT.SauviacL.CottretL.et al (2014). An integrated analysis of plant and bacterial gene expression in symbiotic root nodules using laser-capture microdissection coupled to RNA sequencing.Plant J.77817–837. 10.1111/tpj.12442

  • 51

    RoyS.LiuW.NandetyR. S.CrookA.MysoreK. S.PislariuC. I.et al (2020). Celebrating 20 Years of Genetic Discoveries in Legume Nodulation and Symbiotic Nitrogen Fixation.Plant Cell3215–41. 10.1105/tpc.19.00279

  • 52

    SaitouN.NeiM. (1987). The neighbor-joining method: a new method for reconstructing phylogenetic trees.Mol. Biol. Evol.4406–425.

  • 53

    SerovaT. A.TsyganovaA. V.TsyganovV. E. (2018). Early nodule senescence is activated in symbiotic mutants of pea (Pisum sativum L.) forming ineffective nodules blocked at different nodule developmental stages.Protoplasma2551443–1459. 10.1007/s00709-018-1246-9

  • 54

    TittabutrP.SripakdiS.BoonkerdN.TanthanuchW.MinamisawaK.TeaumroongN. (2015). Possible Role of 1-Aminocyclopropane-1-Carboxylate (ACC) Deaminase Activity of Sinorhizobium sp. BL3 on Symbiosis with Mung Bean and Determinate Nodule Senescence.Microbes Environ.30310–320. 10.1264/jsme2.ME15120

  • 55

    VasseJ.De BillyF.CamutS.TruchetG. (1990). Correlation between ultrastructural differentiation of bacteroids and nitrogen fixation in alfalfa nodules.J. Bacteriol.1724295–4306. 10.1128/jb.172.8.4295-4306.1990

  • 56

    VercruysseM.FauvartM.BeullensS.BraekenK.ClootsL.EngelenK.et al (2011). A comparative transcriptome analysis of Rhizobium etli bacteroids: specific gene expression during symbiotic nongrowth.Mol. Plant Microbe Interact.241553–1561. 10.1094/MPMI-05-11-0140

  • 57

    VoroshilovaV. A.BoestenB.TsyganovV. E.BorisovA. Y.TikhonovichI. A.PrieferU. B. (2001). Effect of mutations in Pisum sativum L. genes blocking different stages of nodule development on the expression of late symbiotic genes in Rhizobium leguminosarum bv. viciae.Mol. Plant Microbe Interact.14471–476. 10.1094/MPMI.2001.14.4.471

  • 58

    WangS.HaoB.LiJ.GuH.PengJ.XieF.et al (2014). Whole-genome sequencing of Mesorhizobium huakuii 7653R provides molecular insights into host specificity and symbiosis island dynamics.BMC Genomics15:440. 10.1186/1471-2164-15-440

  • 59

    YangL.El MsehliS.BenyaminaS.LambertA.HopkinsJ.CazarethJ.et al (2020). Glutathione Deficiency in Sinorhizobium meliloti Does Not Impair Bacteroid Differentiation But Induces Early Senescence in the Interaction With Medicago truncatula.Front. Plant Sci.11:137. 10.3389/fpls.2020.00137

  • 60

    YuanS. L.LiR.ChenH. F.ZhangC. J.ChenL. M.HaoQ. N.et al (2017). RNA-Seq analysis of nodule development at five different developmental stages of soybean (Glycine max) inoculated with Bradyrhizobium japonicum strain 113-2.Sci. Rep.7:42248. 10.1038/srep42248

  • 61

    ZhangL.LiuJ. Y.GuH.DuY.ZuoJ. F.ZhangZ.et al (2018). Bradyrhizobium diazoefficiens USDA 110- Glycine max Interactome Provides Candidate Proteins Associated with Symbiosis.J. Proteome Res.173061–3074. 10.1021/acs.jproteome.8b00209

  • 62

    ZhouS.ZhangC.HuangY.ChenH.YuanS.ZhouX. (2021). Characteristics and Research Progress of Legume Nodule Senescence.Plants10:1103. 10.3390/plants10061103

Summary

Keywords

soybean symbiotic nitrogen fixation, soybean developmental stages, rhizobial gene expression, rhizobial DEGs, nodule development and senescence

Citation

Yuan S, Zhou S, Feng Y, Zhang C, Huang Y, Shan Z, Chen S, Guo W, Yang H, Yang Z, Qiu D, Chen H and Zhou X (2021) Identification of the Important Genes of Bradyrhizobium diazoefficiens 113-2 Involved in Soybean Nodule Development and Senescence. Front. Microbiol. 12:754837. doi: 10.3389/fmicb.2021.754837

Received

07 August 2021

Accepted

04 October 2021

Published

11 November 2021

Volume

12 - 2021

Edited by

Alok Kumar Srivastava, National Bureau of Agriculturally Important Microorganisms (ICAR), India

Reviewed by

Yuan-Ming Zhang, Huazhong Agricultural University, China; Esther Menendez, University of Evora, Portugal

Updates

Copyright

*Correspondence: Haifeng Chen, Xinan Zhou,

†These authors have contributed equally to this work

This article was submitted to Microbial Symbioses, a section of the journal Frontiers in Microbiology

Disclaimer

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.

Outline

Figures

Cite article

Copy to clipboard


Export citation file


Share article

Article metrics