Microbiome analysis revealed distinct microbial communities occupying different sized nodules in field-grown peanut

Legume nodulation is the powerhouse of biological nitrogen fixation (BNF) where host-specific rhizobia dominate the nodule microbiome. However, other rhizobial or non-rhizobial inhabitants can also colonize legume nodules, and it is unclear how these bacteria interact, compete, or combinedly function in the nodule microbiome. Under such context, to test this hypothesis, we conducted 16S-rRNA based nodule microbiome sequencing to characterize microbial communities in two distinct sized nodules from field-grown peanuts inoculated with a commercial inoculum. We found that microbial communities diverged drastically in the two types of peanut nodules (big and small). Core microbial analysis revealed that the big nodules were inhabited by Bradyrhizobium, which dominated composition (>99%) throughout the plant life cycle. Surprisingly, we observed that in addition to Bradyrhizobium, the small nodules harbored a diverse set of bacteria (~31%) that were not present in big nodules. Notably, these initially less dominant bacteria gradually dominated in small nodules during the later plant growth phases, which suggested that native microbial communities competed with the commercial inoculum in the small nodules only. Conversely, negligible or no competition was observed in the big nodules. Based on the prediction of KEGG pathway analysis for N and P cycling genes and the presence of diverse genera in the small nodules, we foresee great potential of future studies of these microbial communities which may be crucial for peanut growth and development and/or protecting host plants from various biotic and abiotic stresses.


Introduction
Bulk and rhizosphere soil contain an array of microbes including bacteria, fungi, algae, protozoa, and archaea (Msimbira and Smith, 2020) with great influence on plant health and nutrient cycling. However, the proper functioning of soil microbial communities depends on many factors including the availability of nutrients, temperature, water, soil structure, pH, and host genotypes and their secreted root exudates (Zgadzaj et al., 2016;Li et al., 2017;Pramanik et al., 2020;Vaidya and Stinchcombe, 2020). An understanding of endophytic or epiphytic soilmicrobe interactions with host roots is of great importance because of their potential to promote plant growth in various ways through enhanced nutrient, mineral and water use efficiency, increased production of phytohormones, induced biotic and abiotic stress tolerance, and minimized toxicity of contaminants (Bacon and White, 2016;Kumar et al., 2017;Harman and Uphoff, 2019;Bhattacharyya and Furtak, 2022).
Studies have shown that nitrogen (N 2 ) fixation has been a great interest for decades because of the advantage of mutualistic soil rhizobia interactions with the host legumes, thereby developing a symbiotic organ "nodule" (Oldroyd, 2013;Liu and Murray, 2016;Zipfel and Oldroyd, 2017;Lace and Ott, 2018;Poole et al., 2018). Rhizobial bacteria housed in this "nodule" fix atmospheric N 2 into ammonia (Martínez-Hidalgo and Hirsch, 2017), a source of N for plant growth and development. The nodule formation and development process in legumes is quite complex and highly restricted due to host specificity (Dazzo and Gardiol, 1984;Menéndez et al., 2019;Walker et al., 2020). Unlike the symbiotic host specificity of the nodulation process, a variety of non-rhizobial endophytic or epiphytic strains can mutually interact with both legumes and non-legume crops (Chen et al., 2003;Berg, 2009;Aserse et al., 2013;Martínez-Hidalgo and Hirsch, 2017;Harman and Uphoff, 2019;Olanrewaju and Babalola, 2019;Fadiji et al., 2020;Swarnalakshmi et al., 2020). Some recent microbiome studies have identified that both rhizobial and non-rhizobial strains can eventually nodulate legumes in which the previously known legume nodulation concept of host specificity/ restriction for rhizobial interactions changed dramatically (Moulin et al., 2001;Sy et al., 2001;Martínez-Hidalgo and Hirsch, 2017;Mayhood and Mirza, 2021;Bender et al., 2022).
Like other legumes, peanuts can form N 2 -fixing root nodules interacting with host-specific rhizobia in cultivated soils. In general, most legume genera (~ 75%) utilize root-hair dependent infection processes (Quilbé et al., 2022); however, the rhizobial infection process in peanut along with some other legumes (such as Sesbania, Aeschynomene, Mimosa, etc.) is unique and occurs through a "crack entry" mechanism (Boogerd and van Rossum, 1997;Held et al., 2010;Sharma et al., 2020), where rhizobia enter through epidermal cell layer of root tissues at the base of emerging lateral roots and proliferation later occurs into the root cortical cells, subsequently developing nodule primordia and nodules (Bogino et al., 2011;Nievas et al., 2012). Previously published research reported that nodule formation in peanut can occur with a broad range of Bradyrhizobium and other rhizobial strains, isolated from different groups of legume plants (Tatiana et al., 2003;Taurian et al., 2006;El-Akhal et al., 2008;Yang and Zhou, 2008;Zaiya Zazou et al., 2018). However, nodulation and N 2 fixation in peanut can vary depending on numerous factors, such as strains and inoculation practices, geographical locations, soil structure, temperature, water availability, soil acidity, salinity, intercropping, crop rotation, fertilization regime, application of herbicide and pesticide uses, etc. (Rowland et al., 2015).
In addition to the above-mentioned factors, it has been reported that nodule size and root zone depth under the soil (particularly during stress condition) have a direct impact on nodule N 2 fixation, effectiveness, and activity (Hardarson et al., 1989;King and Purcell, 2001;Voisin, 2003;Tajima et al., 2007). Tajima et al. (2007) found that N 2 -fixing activity directly correlates with the nodule size. For example, medium size nodules (1.5-2.0 mm in diameter) had the highest N 2fixing activity per unit fresh weight of nodule compared with the smaller (<1.5 mm) or larger (>2.0 mm) size nodules. In general, nodule formation in legumes is variable and their N 2 -fixing activity may change with developmental stages (Tajima et al., 2007). During early growth stages when nodules are smaller in size, typically, the nodule interior color is white and N 2 -fixing activity is low due to the bacterial cell growth and development. However, at the flowering stage when the nodule color is red/pink and medium in size, N 2 -fixing activity is higher due to active enzyme activity and proper bacterial cell growth. At the late stages of plant growth, although the nodule size may be larger, the nodule color changes to greenish and N 2 -fixing activity rapidly decreases. In addition to the nodule size, effectiveness, and N 2 -fixing activity, the lack of resource (C-source) allocation/or sanctioning by plants to the strains that colonize the nodule can also impact nodule size, development, N 2 fixation activity due to bacterial competition, and/or fitness (Westhoek et al., 2021).
In previous studies (Ellman-Stortz, 2021), we noticed that peanut roots in experiments at AgriLife research station (Vernon, Texas, United States) formed two types of nodules: regular size N 2 -fixing nodules (hereafter referred to as big nodules) as well as many small nodules in the same plant root system. Interestingly, the growth of these small nodules was arrested yet the nodule development persisted at all growth stages throughout the plant life cycle. Hence, our main objective and research interest was to understand what microbial communities inhabit these small nodules, which might provide a new direction of research opportunities for peanut nodule function and crop improvement.

Field and experimental conditions
The study was conducted under furrow irrigation at the Texas A&M AgriLife Research and Extension Center at Vernon. The soil is classified as a Miles fine sandy loam (Fine-loamy, mixed, super active, thermic Typic Paleustalfs). Peanuts (Span17) were planted on 14 May 2021 at 16 seed m −1 using a four-row vacuum planter (1-m row spacing). The preceding year's crop was cotton. Seed was treated with a commercial inoculum (Exceed ® Peat for peanut/cowpea/lespedeza/ mung bean; Visjon Biologics) prior to planting. Plots were arranged in a randomized complete block design and replicated four times. Treatments included no cover crop and various cover crop treatments. For this study, peanut plants were collected from the no cover crop control plots only.
Sample harvesting, surface sterilization, and DNA extraction Nodules of two independent plant roots were collected from randomized plots in three developmental stages (R2, R4 and R7; Prasad, 1999). Whole peanut plants with roots and nodules were transported overnight on ice from the field site to the laboratory. Roots and nodules were washed thoroughly with tap water and wiped with a paper towel before collecting two types of nodules (big and small) from roots for surface sterilization. Ten big nodules (~ 1.5-2.0 mm), five from each plant, and 20 small nodules (<1.5 mm), 10 from each plant, were selected and immersed in 70% (v/v) ethanol for 2 min, washed with sterile H 2 O three times, and then soaked in 10% (v/v) commercial bleach for 5 min. After bleach treatment, Frontiers in Microbiology 03 frontiersin.org nodules were carefully washed 5-6 times with sterile H 2 O and rinsed of any trace elements of bleach. A total of 12 big and 12 small nodule samples (3 growth stages x 4 plots) were independently collected for total DNA isolation. Total DNA was extracted from the nodules using DNeasy Plant Pro Kit (Qiagen), according to the manufacturer's instructions. DNA quality and concentration were determined using a NanoDrop 2000 spectrophotometer (Thermo Scientific, Waltham, United States).

Amplification of nifH, nodC, and 16S rRNA genes
To confirm DNA amplifiability, the 16S rRNA gene was amplified using primers Eub338 and Eub518 (Fierer et al., 2005). To check for the presence of symbiotic marker genes (nifH and nodC) in the isolated nodules DNA, primers for nifH (Poly et al., 2001) and nodC genes (Sarita et al., 2005) were used for PCR amplification. The same PCR conditions were used for all three amplicons (16S, nifH, and nodC) which is as follows: 1 cycle of 95°C for 5 min (preheating), 30 cycles of 95°C for 30 s, 55°C for 30 s, and 72°C for 1 min, and a final extension period at 72°C for 5 min. The size of the PCR product was confirmed by electrophoresis on a 2% agarose gel. For PCR reactions, nuclease-free PCR grade water, 2x GoTaq Master mix, and polymerase were purchased from Promega.

Library preparation and sequencing
DNA sequencing was performed at the Texas A&M University Institute for Genome Sciences and Society (TIGSS) using Illumina MiSeq. Amplicon of V3-V4 regions of 16S rRNA gene amplified from 10 ng of DNA and 250 bp paired-end library cleaning, indexing, library quantification, normalization, pooling, and loading for MiSeq sequencing was carried out following the standard library preparation guide (Illumina, United States).
After analyzing raw sequence data, we found more than 5 million sequence reads for both forward and reverse fragments from the 24 nodule microbiome samples with an average count of 111,646 reads per sample generated (Supplementary Table 1). The demultiplexed sequence data resulted an average of 59,560 non-chimeric sequence reads (52.66%; Supplementary Table 2) and inferred a total of 950 amplified sequence variants (ASVs) with an average number of features 40 and mean frequency 59,812 per sample. Since one of the samples (05-PNMB-104R4; Supplementary Table 2) generated a smaller number of non-chimeric sequence reads (6.71%), we further filtered out ASVs based on sampling depth to those that were less than 3,281 feature count/ samples, and the minimum feature frequency was set to 10, leaving behind 492 ASVs from 23 samples (lowest feature count was retained to 35,401), which was analyzed further for community diversity and taxonomic classification.
A taxonomic bar-plot was generated using the filtered table and assigned taxonomy file and visualized through the QIIME 2 view website. The resulting filtered sequences from above were used to construct a phylogenetic tree using "align-to-tree-mafft-fasttree" pipeline from the q2-phylogeny plugin. 3 Subsequently, phylogenetic diversity metrics were constructed with the resulting output "rooted phylogenetic tree. " Alpha and beta diversities were calculated through the q2-diversity plugin 4 using the "core-metrics-phylogenetic" method with a sampling depth (rarefaction) of 1,290. Associations between categorical metadata columns (two types of nodule samples) and alpha diversity data were calculated through "qiime diversity alphagroup-significance" with Kruskal-Wallis pairwise test (Kruskal and Wallis, 1952) for Shannon's diversity index. In addition, we also used PERMANOVA (default) to evaluate if the distances between the two nodule groups (sample compositions) were similar using the "qiime diversity beta-group-significance" plugin (Anderson, 2001) for weighted UniFrac distance. Further, we also conducted Mann-Whitney/Wilcoxon significance tests for alpha and beta diversity in two nodule groups and Kruskal-Wallis pairwise test for beta diversity analysis in three developmental growth stages of small nodules using the "rstatix" package in R. To find similarity within samples or dissimilarity between sample groups, a principal component analysis (PCA) was conducted using the centered log ratio (CLR) transformed microbial compositions. The PERMANOVA (n = 999) was done on Aitchison distance matrices to test significant effect of microbial compositions between two nodule types using the vegan package in R. A compositional heatmap of CLR-transformed microbial abundance (prevalence = 5) across samples was made with R package "microViz. " Frontiers in Microbiology 04 frontiersin.org To identify differentially abundant core microbial communities from two types of nodule samples, we used a statistical power analysis called "ANCOM" (analysis of composition of microbiomes) using ANCOM plugin (Mandal et al., 2015) with "qiime composition ancom" command. For ANCOM analysis, the ASV table that filtered out for low features, frequency, and for chloroplast, mitochondria, Archaea, and Eukaryota was used. Initially, taxa were collapsed using "qiime taxa collapse" command with --p-level 6 to create a genus level table and later, add-pseudocount 1 for log0 (if any), and finally, "qiime composition ancom" command was used to tabulate significant (W-score) groups of bacteria at genus level that were differentially abundant (Estaki et al., 2020).
To predict functional profiles of pathways from 16S rRNA gene sequences of big and small nodule microbiome data (i.e., ASVs and abundance), we analyzed functional metabolic pathways for KO, EC, TIGRFAM, COG, CAZymes, and Pfam using a recently published metagenome analysis tool called "MicFunPred" (Mongad et al., 2021).

Data availability statement
Raw sequencing data were submitted to the NCBI Sequence Read Archive (SRA) under Bio Project accession PRJNA865795 and will be available upon publication. Data files and scripts that were used to generate figures are available at github.com/shak71/PNM.

Peanut forms a distinct pattern of nodulation in the field
Despite being inoculated with host-specific rhizobia, the fieldgrown peanuts formed both regular (referred to as big) sized ( Figure 1A) and small nodules ( Figure 1B; Please see M&M section for nodule size). These small nodules were present at all observed growth stages throughout the plant life cycle. While the larger nodules had a red interior color suggesting active N 2 fixation, the small nodules showed variable color (white, light pink, and gray or green) suggesting limited N 2 fixation activity (data not shown). However, N 2 fixation genes were found in both types of nodules. The nifH gene was found in all big and small nodules, except one small nodule while the nodC gene was not amplified from few small nodule samples (Supplementary Table 3).

Diversity analysis of bacterial communities in peanut nodules
To assess and investigate the bacterial community richness and diversity in two nodule types, the alpha diversity measures for observed ASVs and Shannon Diversity Index were calculated in big and small nodules of peanut roots. The number of observed ASVs differed significantly (p = <0.0001) between the two nodule types (Figure 2A). Observed ASVs in big nodule samples ranged from 5 to 10 while in small nodule samples ranged from 11 to 235. A Mann-Whitney/Wilcoxon significance test analysis for Shannon Diversity index ( Figure 2B) was used to find the community richness which indicated a significant difference and variation between the big and small nodule microbiomes (p = 0.00029) with the bacterial communities in the small nodules being more diverse than those in the big nodules.
To analyze the bacterial community compositions in two types of nodules as well as growth stages for small nodules, the beta diversity index was calculated based on Weighted UniFrac distances. A Mann-Whitney/Wilcoxon significance test analysis was used among all samples grouped as big and small for Weighted UniFrac distance for the significance test and the results clearly showed that bacterial compositions differ between big and small nodule types (p = < 0.0001; Figure 2C). Since bacterial compositions differ within the small nodule category, we further explored the possibility within three developmental growth stages of small nodules and indeed, we found that bacterial compositions differ significantly (Kruskal-Wallis test, p = <0.0001) among all three growth stages of small nodules ( Figure 2D). Further compositional aspects of centered log ratio (CLR) transformed PCA were analyzed for the bacterial communities in the big and small nodules as well as three different growth stages (R2, R4, and R7) of the small nodules. Based on the PERMANOVA test, the PCA results indicated that there were significant differences in bacterial communities between big and small nodules (p = 0.001); however, the bacterial communities among three growth stages of small nodules were also slightly different (p = 0.061; Figure 3A). Within small nodules, bacterial communities tended to group separately in these three growth stages except one or two data points in R4 and R7 stages that were closer to other growth stages. In contrast, bacterial communities in the big nodules seemed to cluster together and vary less across growth stages. Notably, the small nodule communities were initially (R2) somewhat similar to the big nodules but became more distinct at later growth stages (R4 & R7). At the genus level, based on top taxa loading on the PCA plot as well as direct visualization of sample compositions on the iris plot ( Figures 3A,B), big nodule samples contained mostly Bradyrhizobium, while other taxa were often more abundant in small nodule samples. Collectively, the data clearly indicated that the small nodules had more diverse and temporally variable bacterial communities than the big nodules.

Bacterial community composition in big and small nodules of peanut roots
All the filtered ASV's from the big and small nodules were classified as bacteria and included 15 phyla, 25 classes, 57 orders, 87 families, and 151 genera (Figure 4; Supplementary Table 4). Both nodule types had large relative populations of Proteobacteria ( Figure 4A), with the phylum representing ~99.96% of the bacterial community in the big nodules and ~ 68.67% in the small nodules combinedly from three growth stages. Thus, the big nodules harbored a tiny fraction of native community microbes while the small nodules were inhabited by a large portion (~31.32%) of diverse microbes. Further analysis revealed that only Actinobacteriota (~0.03%) co-inhabited with Proteobacteria in the big nodules; however, 15 Frontiers in Microbiology 05 frontiersin.org different phyla including Actinobacteriota were detected in the small nodules ( Figure 4B). Interestingly, the microbiome data showed that these diverse bacterial populations in the small nodules dominated over time, becoming more abundant at later growth stages -the most noticeable being members of the Actinobacteriota, Bacteroidota, Chloroflexi, and Patescibacteria ( Figure 4A). Surprisingly, the phyla Patescibacteria and Chloroflexi were not enriched at all in the early R2 growth stage. In addition, we found that Proteobacteria abundance and occupancy in the big nodules was steady throughout the growth stages. However, Proteobacteria were comparatively lower while the diverse microbial populations appeared higher in the small nodules irrespective of growth stages. At the phylum level, the microbial abundance and distribution frequency was analyzed since the most dominant and abundant phyla were Proteobacteria and Actinobacteriota compared to other taxa. The data revealed that microbial abundance and frequency varied from low to high abundance in the big and small nodules for these two phyla with more distinct bacterial populations present in small nodules. To further explore the variation of microbiota within these three growth stages of the small nodules, 118 taxa at the genus level were extracted based on relative abundance without any uncultured taxa ( Figure 4C; Supplementary Table 5). The Venn diagram showed the taxa (relative abundance >0) distribution in R2, R4 and R7 growth stages of small nodules, and indicated there was diverse microbial dominance over the growth cycle of the small nodules. Interestingly, we also found some unique as well as overlapped microbiota between these three growth stages with the greatest overlap between growth stages R4 and R7.
To identify the dominant and abundant taxa in the big and small nodules, we first removed all the unclassified species-level data, the rare and uncultured taxa, the taxa with a minimum prevalence of 5, and then a CLR-transformed compositional heatmap was made at the genus level ( Figure 5A). The small nodule samples from the three growth stages had a diverse bacterial community while the big nodule samples were dominated by a single genus Bradyrhizobium spp. The core bacterial genera in the big and small nodules were identified using ANCOM differential test analysis at the genus level ( Figure 5B). Based on 'W' score > 100, 7 genera were identified: Bradyrhizobium, Streptomyces, Niastella, Allorhizobium-Neorhizobium-Pararhizobium-Rhizobium, Sphingomonas, Novosphingobium, and Chitinophaga. These seven core genera were also found in the compositional heatmap, suggesting that these are the dominant and highly abundant core genera in big and small nodules. Of these seven abundant and dominant core genera, based on W score, only Bradyrhizobium was highly abundant in the big nodules, while at least 6 genera including Streptomyces, Allorhizobium-Neorhizobium-Pararhizobium-Rhizobium, Niastella, Novosphingobium, Sphingomonas, and Chitinophaga were primarily abundant in the small nodules.

Prediction of functional pathway genes of bacterial communities
Predicted functional profiles of pathway genes (based on pathway abundance) of bacterial communities from big and small nodules of peanut roots were generated using 'MicFunPred' (Mongad et al., 2021). Utilizing 16S rRNA gene sequences and normalized ASV's abundance from three developmental stages of big and small nodules, metagenomes in terms of KEGG Orthology (KO) functional categories (based on log10 transformed pathways abundance) were used to investigate the enriched function for pathway genes related to nitrogen (N) and phosphorus (P) cycling processes (Figures 6, 7; Supplementary Table 6). Based on search criteria as "nitrogen metabolism" at C-level hierarchy (Supplementary Table 6), a total of 32 genes from the big nodules and 36 genes from the small nodules related to the N cycling process were enriched. A heatmap was made based on all the 32 genes that were common in both big and small nodule developmental stages ( Figure 6). The relative abundance of most biological N 2 fixation (nifH, nifD, and nifK), assimilatory (nirA, nasA, and nasB) and dissimilatory (nirB, napA, and napB) nitrate reduction, and denitrification (norB and norC)-related genes were higher in all three growth stages of big nodules compared to small nodules at all growth stages ( Figure 6). However, the predicted pathway data showed that the relative abundance of some genes or different subunits related to the N cycling process was higher in at Field-grown peanut root-nodules used as a source of DNA for 16S-based microbiome sequencing and analysis. (A,B) represent big nodules and small nodules, respectively. The images were taken with the same magnification (scale 50 mm).  Figure 6). Interestingly, a few genes related to N cycling processes that were uniquely enriched only in the small nodules were observed. These unique nitrogen-related pathway genes were Glutamate synthase (ferredoxin; GLU, gltS), ferredoxin-nitrate reductase (narB), nitrous-oxide reductase (nosZ), and hydroxylamine reductase (hcp; Supplementary Table 6). Like the N cycling processes, KO pathway genes related to the P cycling processes were analyzed (Figure 7). Based on search criteria "phosphate metabolism/pathway" at C-level hierarchy (Supplementary Table 6), a total of 48 genes from three growth stages of big nodules and 58 genes from small nodules related to P cycling processes were enriched. A heatmap was made based on all the 48 genes that were common in both big and small nodule developmental stages (Figure 7). Among 48 phosphate-related genes, almost 50% (26) genes showed higher relative abundance in big nodules growth stages compared to small nodules. However, the relative abundance of 22 Frontiers in Microbiology 07 frontiersin.org genes related to P cycling processes was also higher in all three stages of small nodules category, particularly 3-phytase and 4-phytase/acid phosphatase, myo-inositol 2-dehydrogenase (iolG), fructose-1,6bisphosphatase II (glpX), deoxyribose-phosphate aldolase (deoC), ribose 5-phosphate isomerase B (rpiB), myo-inositol-1-phosphate synthase (INO1), 6-phosphofructokinase 1 (pfkA), transaldolase (talA, talB), and glucose-6-phosphate isomerase (GPI; Figure 7). The unique phosphaterelated genes enriched in small nodules were inositol oxygenase (MIOX), 1-phosphatidylinositol phosphodiesterase (plc), 2-keto-myoinositol isomerase (ioll), diphosphate-dependent phosphofructokinase (pfp), 6-phospho-3-hexuloisomerase (hxlB), 3-hexulose-6-phosphate synthase (hxlA), dehydrogluconokinase (kguK), glucose-6-phosphate isomerase-archaeal (pgi1), and CDP-diacylglycerol--inositol 3-phosphatidyltransferase (CDIPT), and fructose-1,6-bisphosphatase III (fbp3). Besides N & P cycling genes, other pathway genes were analyzed that were unique only in the small nodule growth stages. All the duplicates for KO identifications from both big and small nodules were removed, and then sorted out for unique pathway genes in small nodule category only, of which 383 genes were observed that were unique and enriched (log10 abundance >4.0 in all three stages; Supplementary Table 7). Among the observed 383 pathway genes, a heatmap was made for 35 genes that were common among three stages of pre-selected top 50 pathway genes (Figure 8). All 35 pathway genes were highly enriched in growth stage R7 compared to R2 and R4, suggesting that these genes might play some role in the late stage of the small nodules. The most noticeable and highly enriched genes were IMP dehydrogenase (IMPDH, guaB), NADH-quinone oxidoreductase subunit J, K, N (nuoJ, nuoK, and nuoN), and succinate dehydrogenase / fumarate reductase (sdhC, frdC), and phosphate regulon sensor histidine kinase PhoR (phoR). It is quite interesting to observe that some key metabolic categories other than N and P were specifically increased in the small nodules (Supplementary Table 7) and this warrants future studies for further characterization.

Discussion
In this study, the nodule microbiome from two distinct types of nodules formed on the same field-grown peanut roots was characterized using a 16S rRNA amplicon strategy. Interestingly, we noticed that the occurrence and abundance of small nodule phenotype was consistent across plant growth stages, as observed in previous studies (Ellman-Stortz, 2021). It is well known that the lack of some host genetic components tightly regulates bacterial infection processes and may control this type of phenotype (Hossain et al., 2012). However, it was observed that the field-grown peanut with the commercial seed coating inoculum formed small nodules and persisted at all growth stages throughout the plant life cycle alongside the regular N 2 -fixing nodulation. There is evidence that small nodules may result from resource sanctioning by the plant in response to specific nodule occupants (Westhoek et al., 2021); however, the microbiome composition of these smaller nodules has not been extensively investigated, especially for peanuts. This was the impetus to understand and classify the inhabitants of these two distinct nodule types.
It was observed that nifH and nodC (symbiotic marker genes) were generally present in both the big and the small nodules, suggesting that bacterial communities in these nodules might have the ability to fix N 2 . We did not directly measure the gene expression or N 2 -fixing activity in these nodules, but the examination of several nodules indicated a light pink/green color in small nodules and pink color in big nodules (data not shown). Also, previous reports suggest that small nodules (<1.0 mm) with light red/pink color might have less N 2 -fixing activity than medium sized (1.5-2.0 mm) nodules (Tajima et al., 2007). The non-detection of nifH or nodC gene in a few small nodule samples might be due to the low copies of the gene in the genome or the possibility of missing genes from the bacterial genomes. The latter is less likely because the same samples in other replications and/or growth stages showed the presence of the gene (Supplementary Table 3).
After utilizing filtration criteria, our 16S rRNA analyses from nodules revealed only 492 ASV's. The big nodules had less observed ASV's with a higher frequency while small nodules possessed more ASV's compared to big nodules with comparatively less frequency, which further indicates a tiny fraction of community diversity/ Principal Component Analysis (PCA) and sample compositions of bacterial communities in big and small nodules of peanut roots. (A) Centered log-ratio transformed (CLR) an unconstrained genus level PCA analysis in peanuts big (11) and small (12) nodule samples and three developmental growth stages (R2, R4, R7). The top five taxa (labeled with their annotations) were indicated by bold black lines. (B) Iris plot showing the sample compositions of microbiota at the genus level in big and small nodule growth stages (R2, R4, and R7). The big (B) and small (S) nodule samples on the iris plot were automatically arranged by their rotational position around the center/origin of the PCA plot. Both PCA and iris plots were performed with R package "microViz." Frontiers in Microbiology 08 frontiersin.org occupants in big nodules. Since most of the ASV's came from small nodules, initial analysis of these ASV's clearly indicates that small nodule occupants were more diverse than big nodules, which might be the reason for the bimodal distribution patterns. Further analysis of this data through taxonomic community profiling found more diverse bacterial communities occupied in the small nodules compared to the big nodules, suggesting that these diverse microbial communities in the small nodule might play some role(s) for microbemicrobe and/or microbe-plant interactions for the benefit of plant Relative abundance of bacterial communities at phylum level in (A) big and small nodules of peanut roots, (B) distribution and relative abundance of microbial communities of two major phyla Actinobacteriota and Proteobacteria, and (C) a Venn diagram showing the taxa at genus level from three developmental growth stages of small nodule type. For Venn diagram, taxa were selected based on relative abundance >0 for three growth stages of small nodules (R2, R4 and R7), and was made using Venny 2.1.0 (https:// bioinfogp.cnb.csic.es/tools/venny/).

FIGURE 5
Microbial abundance and the differentially abundant core bacteria in big and small nodules of peanut roots. (A) Heatmap showed the centered log-ratio transformed (CLR) microbial abundance at the genus level. The x and y-axis indicate sample name and the microbial taxa, respectively. Peanut nodule microbiome big (PNMB) and small (PNMS) denoted on the x-axis. 101, 202, 301, and 401 represent biological replications, and R2, R4, and R7 represent three developmental growth stages. (B) The plot showed the differentially abundant core microbiota at the genus level in big and small nodules. The differentially abundant core microbiota was selected based on ANCOM analysis as stated in the Materials and Methods section. The W value represents the number of times the null hypothesis (the average abundance of a given genus in a group is equal to that in the other group) was rejected for a given genus. The x-axis value represents the CLR-transformed mean difference in abundance of a given genus between the big and small nodule groups. A positive x-axis indicates a genus is abundant in a small nodule group and a negative x-axis value indicates a genus is abundant in a big nodule group. Only genera with W value >20 were plotted (color coded).
Frontiers in Microbiology 09 frontiersin.org growth promotion, microbe/plant adaptation to harsh environmental stress, and/or protecting from disease or environmental pollution. It is noteworthy to mention that peanut seeds used for the experiment were coated with commercial inoculum before planting. Since rhizobia inoculants often fail to compete for nodule occupancy against the native soil rhizobia due to the strain's variation, abundance, and composition as well as many environmental factors (Ji et al., 2017;Mendoza-Suárez et al., 2021), we had the opportunity to look at the possibility of rhizobia competition and occupancy in our data set. In the context of competition and occupancy prospective, we did not see any major evidence for native microbial community members displacing the commercial inoculum in the big nodules, suggesting that commercial inoculum might be active in big nodules throughout the plant growth cycle. However, it was quite surprising to observe that a diverse set of native microbial community members occupied the small nodules. These native communities gradually dominated over time in the small nodule compartments but did not have any clear indication about why native microbial communities colonized small nodules from the early stage and kept competing throughout the peanut life cycle. However, our data from microbial composition analysis as well as core dominant microbial community profiling of two types of nodule categories provide some intriguing insights into this. Based on both microbial composition analysis and core dominant microbial genera observed, only Bradyrhizobium was abundant in the big nodules while Streptomyces were primarily abundant & dominant in the small nodules, along with several core genera including -Neorhizobium-Pararhizobium-Rhizobium, Niastella, Novosphingobium, Sphingomonas, Chitinophaga, and Caulobacter, that solely occupied the small nodules.

Allorhizobium
Though we could not directly correlate biological and physiological functions of these core bacteria that occupied the small nodules either as beneficial or commensal and/or pathogenic (Brader et al., 2017), some previous reports highlighted the beneficial role and importance of these endophytic bacterial taxa (Aserse et al., 2013;Khan et al., 2014;Rodriguez-Conde et al., 2016;McKee et al., 2019;Olanrewaju and Babalola, 2019;Asaf et al., 2020;Segura et al., 2021). For example, the genus Niastella, from the Chitinophagales order under phylum Bacteroidota, can degrade a variety of carbohydratebased biomass, releasing sugars as nutrients for their own growth as well as for other microbial communities (McKee et al., 2019). Recent studies found that bacterial species from Novosphingobium, Sphingomonas, and Niastella genera play significant roles: from remediation of environmental contamination to producing highly beneficial phytohormones (gibberellins and indole acetic acid, sphingan, and gellan gum), improving plant-growth attributes (i.e., shoot length, chlorophyll contents, and shoot and root dry weight) protecting plants from stress conditions such as drought, salinity, and heavy metals in agriculture soils, co-metabolizing nutrients existing in root exudates, and producing the N-acyl-homoserine lactone quorum-sensing (QS) signals (Gan et al., 2009;Khan et al., 2014;Rodriguez-Conde et al., 2016;Asaf et al., 2017Asaf et al., , 2020Segura et al., 2021). Heatmap of the predicted KEGG Orthology (KO) pathway genes of bacterial communities from three developmental stages of big and small nodules in peanut roots related to N cycling processes. The scale bar indicates the color saturation gradient based on log10 transformed values of pathway abundances enriched from bacterial taxa in big and small nodules. R2, R4, and R7 represent three developmental growth stages. The partial data for genes related to N cycling processes used to plot the heatmap and the full set of enriched pathways data is listed in Supplementary Table 6.
Frontiers in Microbiology 10 frontiersin.org The genus Streptomyces under phylum Actinobacteriota was almost exclusively detected in the small nodules. Most species in this group are widely distributed as soil microbes, are efficient as rhizosphere, rhizoplane, and endophytic colonizers of host plants, and are a useful source for producing bioactive compounds, antibiotics, and extracellular enzymes. Thus, this group of bacteria functions for plant growth promotions as a biocontrol agent as well as a biofertilizer. As biocontrol agents, they synthesize plant growth regulators, siderophore production, antibiotic production, and volatile compound secretion. As biofertilizer, they directly promote plant growth promotion by the production of phytohormones (such as auxins, cytokinins, and gibberellins), the scavenging of ferric iron from the environment by siderophores, N 2 fixation, and the suppression of stress response in plants by production of 1-aminocyclopropane-1carboxylate (ACC) deaminase activity (for details please see the review Olanrewaju et al., 2019).
Although there are several tools available to predict functional profiles of 16S rRNA gene sequence data (DeSantis et al., 2006;Pruesse et al., 2007;Douglas et al., 2020;Narayan et al., 2020;Wemheuer et al., 2020), a recently developed bioinformatic tool "MicFunPred" (Mongad et al., 2021) was utilized because of its advantages using a novel approach to predict and understand the functional profiles of the metagenome based on a set of core genes only, minimizing false-positive predictions in the microbial communities of big and small nodules of the peanut roots. Here, the KO functional category for genes related to N and P cycling processes were reported. Relative abundance of many genes were observed related to the N and P cycling processes, which were higher in abundance in big nodules compared to small nodules. However, it was quite interesting that a few genes/subunits showed clearly higher abundance in the small nodules compared to the big nodules. For the N cycling process in small nodule microbiota, nitrate reductase/nitrate FIGURE 7 Heatmap of the predicted KEGG Orthology (KO) pathway genes of bacterial communities from three developmental stages of big and small nodules of peanut roots related to P cycling processes. The scale bar indicates the color saturation gradient based on log10 transformed values of pathway abundances enriched from bacterial taxa in peanut big and small nodules. R2, R4, and R7 represent three developmental growth stages. The partial data for genes related to P cycling processes used to plot the heatmap and the full set of enriched pathways data is listed in Supplementary Table 6.
Frontiers in Microbiology 11 frontiersin.org oxidoreductase alpha, beta, and gamma subunits (narG, narH, and narl), carbamate kinase (arcC), and glutamate dehydrogenase (NADP + ; gdhA; GLUD1_2) were high in abundance. Nitrate reductases reduce nitrate by microbial dissimilatory processes for denitrification or dissimilatory reduction of nitrate to ammonium. The presence of these subunit/genes with higher abundance in the small nodule microbiota suggest they might participate in the N 2 -fixation process. Another example is that glutamate dehydrogenase (NADP + ; gdhA; GLUD1_2) was higher in the small nodule microbiota, while glutamate synthase (NADPH; gltB, and gltD) was higher in the big nodule microbiota. It was previously reported that bacteria or yeasts express distinct GDH isoenzymes for metabolic or biosynthetic processes. Based on nutrient type and availability, organisms that live in an amino acid rich environment use NAD + − dependent GDH for their catabolic needs. In contrast, organisms utilizing inorganic N 2 (nitrates or ammonia) use the NADP + − specific GDH for their synthetic needs (Plaitakis et al., 2017) indicating that microbial communities in big and small nodules might use different GDH subunits for their catabolic or synthetic needs. Like the N cycling process genes, also it was observed that most of the P cycling process genes were higher in abundance in the big nodules microbiota; however, 3-phytase (E3.1.3.8; KO 1083), 4-phytase/ acid phosphatase (appA; KO1093) and phosphogluconate dehydrogenase (edd; KO1690), myo-inositol 2-dehydrogenase (iolG), 6-phosphofructokinase 1 (pfkA), and glucose-6-phosphate isomerase (GPI) were higher in abundance in the small nodules microbiota. There might be some correlation between the presence of microbial communities and the abundance of phytases in the small nodules. Indeed, a recent report suggested that Chitinophagales, Cytophagales, Sphingomonadales, and Xanthomonadales were the order-level microbial indicators that positively correlated with the soil available phosphorus in agricultural soils (Wu et al., 2022). Most likely, this community of microbes capable of producing more phytases for hydrolyzing myo-inositol phosphates (phytates), which are abundant in many soil types, makes up a large portion of bioavailable soil P to improve plant growth and development (Balaban et al., 2017). Besides N-and P-related pathway genes, some other category of pathway genes were observed that specifically were enriched in the small nodules particularly in the late growth stages, suggesting that these pathway genes might play some role in the small nodules as well. The main objective of this study was the community profiling of two distinct field-grown nodule types in peanuts. Though the data could not fully support experimentally for biological or physiological functions for the predicted genes/subunits stated above, future studies will be required to answer this question.

Conclusion
This research clearly observed that microbial inhabitants were quite different between the two types of nodules that formed on the same peanut roots under field conditions. The microbial community in the big nodules were predominantly colonized by the commercial rhizobial inoculum for the formation of N 2 -fixing nodules to fix Heatmap of top 35 predicted KEGG Orthology (KO) pathway genes of bacterial communities that were unique to three developmental stages of small nodules of peanut roots. The scale bar indicates the color saturation gradient based on log10 transformed values of pathway abundances enriched from bacterial taxa in small nodules growth stages. R2, R4, and R7 represent three developmental growth stages.
Frontiers in Microbiology 12 frontiersin.org atmospheric N 2 . In contrast, diverse native microbial communities co-inhabited with the commercial inoculum in the small nodules and dominated gradually during the plant growth and development. Based on our data analysis, the prediction of functional pathway of genes related to the N and P cycling processes as well as with the previous published reports indicated that the diverse microbial inhabitants in the small nodules might play some role in many ways (such as production of phytohormones, bioactive compounds, antibiotics, N and P availability, stress tolerance, phytoremediation for toxic chemicals etc.) for better peanut nodulation, growth, and development. Future studies with our recently isolated endophytes (un-published) from the same peanut nodules used for this microbiome study will certainly verify some of the above-mentioned functional roles and will provide clear data and recommendations for peanut improvement.

Data availability statement
The data presented in the study are deposited in the https://www. ncbi.nlm.nih.gov/ repository, accession number PRJNA865795.

Author contributions
MH designed and planned experimental procedures, collected samples for microbiome sequencing, performed all experimental works, bioinformatics, and statistical analysis, generated figures, and wrote the manuscript. PD designed the experimental field plot, collected and shipped plants to the laboratory from the field site, and revised the manuscript. TG provided advice for experimental design and plan, and edited and revised the manuscript. All authors contributed to the article and approved the submitted version.