Genome-Resolved Metagenomics of Nitrogen Transformations in the Switchgrass Rhizosphere Microbiome on Marginal Lands

: Switchgrass ( Panicum virgatum L.) remains the preeminent American perennial (C4) bioenergy crop for cellulosic ethanol, that could help displace over a quarter of the US current petroleum consumption. Intriguingly, there is often little response to nitrogen fertilizer once stands are established. The rhizosphere microbiome plays a critical role in nitrogen cycling and overall plant nutrient uptake. We used high-throughput metagenomic sequencing to characterize the switchgrass rhizosphere microbial community before and after a nitrogen fertilization event for established stands on marginal land. We examined community structure and bulk metabolic potential, and resolved 29 individual bacteria genomes via metagenomic de novo assembly. Community structure and diversity were not signiﬁcantly different before and after fertilization; however, the bulk metabolic potential of carbohydrate-active enzymes was depleted after fertilization. We resolved 29 metagenomic assembled genomes, including some from the ‘most wanted’ soil taxa such as Verrucomicrobia, Candidate phyla UBA10199, Acidobacteria (rare subgroup 23), Dormibacterota, and the very rare Candidatus Eisenbacteria. The Dormibacterota (formally candidate division AD3) we identiﬁed have the potential for autotrophic CO utilization, which may impact carbon partitioning and storage. Our study also suggests that the rhizosphere microbiome may be involved in providing associative nitrogen ﬁxation (ANF) via the novel diazotroph Janthinobacterium to switchgrass.


Introduction
"Plants wear their guts on the outside" wrote Janzen (1985) [1], since the rhizosphere of terrestrial plants-the~millimeter interface between plant roots and surrounding soilplays critical roles in nutrient uptake, absorption, and degradation, via the diverse microbes it contains [2][3][4][5]. The rhizosphere connects plants to ecosystem processes including the cycling and sequestration of water, nitrogen (N), carbon, and other nutrients [6]. The rhizosphere represents one of the most dynamic and diverse interfaces on the planet, containing up to 10 11 microbial cells per gram root, potentially representing over 30,000 bacterial Switchgrass (Panicum virgatum L.) is the principal United States bioenergy model C4 perennial crop for use in cellulosic ethanol production, biogas, and combustion [32], which could displace up to 30% of the current petroleum consumption [33]. Its high biomass productivity in low-nutrient soils common to marginal lands is key because growth on marginal land avoids competition with food crops on arable lands [34][35][36][37]. Switchgrass and other cellulosic ethanol sources could displace up to 30% of the current petroleum consumption [33].
Contrary to many cropping systems which show N limitation, N fertilization of switchgrass often has limited to no impact on productivity [36,38], and there is a resulting gap in the N budget [36]. This suggests that the rhizosphere microbiome could contribute to the N requirement [39] through free-living diazotrophic bacteria using associative nitrogen fixation (ANF) [25,40,41].
Studies to date have investigated the effect of N fertilizer on the switchgrass rhizosphere microbiome using 16S rRNA amplicon sequencing [39,40,42]. 16S rRNA amplicon studies typically fail to resolve novel 'candidate' phyla, representing vast portions of the tree of life [43,44]. These 16S rRNA amplicon approaches fail to measure metabolic capabilities directly, functional potential, gene-gene assortment on operons via long-range sequence contiguity, gene transfer via horizontal gene transfer (HGT), or the provision of MAGs. Shotgun metagenomics thus provides a powerful approach to characterize both diversity and functionality simultaneously [45].
Here, we use high-throughput metagenomic shotgun sequencing of the switchgrass rhizosphere on marginal lands to resolve its taxonomic composition and functional metabolic potential, and resolve community microbial genomes (MAGs). We further compare plots pre-and post-fertilization to determine how the microbiome responds to N addition and give insight into its responsiveness over a two-week interval. We first assess metagenome quality and overall diversity across our study, then we describe overall metabolic potential and shifts between timepoints, delving particularly into nitrogen cycling (e.g., ANF) and carbohydrate metabolism (i.e., CAZy pathways). Finally, we describe MAGs of abundant bacteria and use these to elucidate their potential roles in nitrogen cycling and coupling the N and C cycles in the rhizosphere.

Study Site Description and Management
Switchgrass rhizosphere soil was sampled (April to May 2016) from the Lux Arbor reserve (42.48 N-85.44 W) as part of the Great Lake Bioenergy Research Center's (GLBRC) marginal land experiments. The Lux Arbor reserve marginal land soil is a sandy clay loam, mesic Hapludalf (Crum and Collins, 1995). The mean annual temperature (MAT) was 10.1 • C, with a mean annual precipitation of 1.005 m. Four blocks of switchgrass rhizosphere soil from the Cave-in-rock variety were sampled in a randomized block design, after a two-week pulse treatment of N amendment. The assigned blocks were 64 ft × 40 ft in size. Amended soil was approximately half a plot, with the other half receiving no amendment. The soil was classified as a Kalamazoo loam with a 2-6% slope. We sampled switchgrass rhizosphere within the G5 blocks, on the block half that was historically amended with N fertilizer, but not at the time of sampling called here pre-fertilized (i.e., I1-I4, prefertilized-initial samples). We sampled the same blocks two weeks post-nitrogen fertilizer amendment, which here is called post-fertilized (i.e., P1-P4, post-fertilized-post samples 1 to 4) ( Figure S1). Blocking diagram and map of the location is provided ( Figure S1).
We collected soil cores (2 cm diameter × 15 cm deep) from near the centers of each halfblock,~50 g of soil per core, one sample taken per half-block. Cores from each half-block on each date were sieved through a 4 mm sieve to remove rocks and large roots, and then flash-frozen using liquid nitrogen until processed. The rhizosphere is typically considered to be the zone of soil that is influenced by roots, so, given that these soil cores were sampled at the base of a well-established perennial plant and that the cores themselves contained a sizeable amount of root material, we considered these samples to be indicative of the switchgrass rhizosphere. The "tightly bound" rhizosphere and rhizoplane are typically defined operationally by the particles of soil that are stuck to a root after shaking and the surface of a root after all soil particles are removed, respectively. We did not seek to partition these compartments, but noted that the tightly bound and rhizoplane microbes were captured in our samples.

DNA Extraction and Sequencing
Total DNA was extracted from~2 g switchgrass rhizosphere from field flash-frozen samples using the MoBio PowerSoil DNA (Carlsbad, CA, USA), according to the manufacturer's instructions. Samples were quantified using the Qubit Fluorometer 2.0 (Invitrogen, Carlsbad, CA, USA) and quality-checked using a Nanodrop-1000 (Thermo Fisher, Waltham, MA, USA). Michigan State University Research Technology Support Facility (RTSF) sequencing core completed Illumina library preparation, library quantification, and sequencing on HiSeq 4000 150 bp paired-end read format. We analyzed eight metagenomes from switchgrass rhizosphere within marginal lands of southern Michigan Lux Arbor Reserve, comprising four samples from each pre-and post-fertilizer application.
DESeq2 R package (version 1.18.1) was used to obtain differential statistics on taxonomic composition and functional annotations from predicted ORFs from contigs using COG, CAZy, and KO annotation abundances [52]. The contig ORF counts' abundances within DESeq2 followed a paired analysis, which blocked by sampling plot pre-/postfertilization, then normalized with variance stabilizing normalization.
We pooled contigs from pre-/post-fertilization, then used the differential read abundances across all samples to obtain MAGs. Contigs were binned using Concoct [53], Maxbin2 [54], and Metabat2 [55], including the refinement program within metawrap [56]. Metawrap was used for refinement of MAGs, blobiology prediction, and quantification of bin (MAG) abundance via quant_bin module [56]. CheckM was used to evaluate completeness, contamination, redundancy, and genome properties of the MAGs [57]. All MAG qualities were reported according to the MIMAG standards [8]. We tested a variety of methods to resolve the taxonomy of the MAGs, including metawrap's classifier, classify genomes (https://github.com/AlessioMilanese/classify-genomes, (accessed on 15 April 2023)), which uses the metagenomic operational taxonomic units (mOTU) v2 taxonomy, JSpeciesWS Tetra Correlation Search [58], GTDB-Tk [59], and blastp of ribosomal protein S9 gene (obtained from Prokka annotation). Only GTDB-Tk provided the taxonomic identifications that were supported by ribosomal protein S9 gene blastp results. GTDB-Tk provided all further taxonomy for the MAGs downstream.

Read-Based mOTU Picking and Statistical Analysis
mOTUs analysis of the quality-filtered and decontaminated reads (not contigs) used in the de novo metagenomic assembly were used as input mOTUs v2 [45], then parsed and further analyzed with the phyloseq (version 1.22.3) R package [60]. Alpha diversity measurements with statistical testing for mOTUs (including t-test, Wilcox, Kruskal, and ANOVA) were completed in the phyloseq (version 1.22.3) R package. Beta-diversity metrics were obtained for mOTUs using the UniFrac (weighted/unweighted) and Brays-Curtis distances in phyloseq (version 1.22.3) R package. DESeq2 R package (version 1.18.1) was used within phyloseq (version 1.22.3) R package for mOTU differential statistics using negative binomial distribution corrected with variance-stabilizing normalization.

Data and Analysis Code Availability
Raw sequence data, assembled contigs, and supplemental data are all available on NCBI Bioproject/BioSample/SRA as PRJNA954399-Lux Arbor Metagenome, and are also available at Open Science Framework, with code, at https://osf.io/mzrvj/ (accessed on 15 April 2023). All code for this study is available on www.github.com/friesenlab/MMPRNT_ panicum_metagenome_mags/ (accessed on 15 April 2023).

Assessment of Assembly and Metagenomic Assembled Genomes within Lux Arbor
The raw data represent 5.37 billion Illumina reads, with 805 Gbp in 490 gigabytes of compressed data, with~100 gigabytes of uncompressed data per sample (Table 1). On average, 34.5% of the data were removed per sample due to quality, length, or adapters, or because they were phiX174 bacteriophage Illumina spike-in DNA library (Table 1). Upon metagenomic de novo assembly with MEGAHIT, each sample averaged 4.6 million total contigs (>200 bp) which contained, on average, 3.2 Gbp, with an average N50 of 737 bp (Table 1). However, the best assessment of a soil and rhizosphere metagenome de novo assembly is the number of contigs which are >1 kbp and are longer than 5 kbp [9,10]. On average,~60,000 contigs per sample were contained on contigs > 1 kbp, with an average of 2.21 Gbp assembly size and an average N50 value of 1982 bp (Table 1). Across all samples, we obtained 190,172 contigs > 5 kbp, with an average of 23,771 contigs of 5 kbp per sample on a 1.19 Gbp assembly size with an average N50 value of 9254 bp (Table 1). Of those 190,172 > 5 kbp contigs, 44,171 of them were >10 kbp in length, and 237 were >100 Kbps. Max contig length was 697,599 across all samples.

Microbiome Diversity and Composition of Lux Arbor Switchgrass Rhizosphere
Only 571 mOTUs were identified across all 8 samples combined. Comparing alpha diversity between pre-and post-fertilization samples showed no statistically significant difference between diversity metrics, which included observed mOTUs, ACE richness, Shannon diversity, or Simpson evenness ( Figure 1A, Table S1). Qualitatively, the variances were more substantial and more variable for alpha diversity metrics within post-fertilized plots than pre-fertilized ( Figure 1A). There was no statistically significant difference between alpha diversity metrics ( Figure 1A) due to high variance observed post-fertilization. Using UniFrac (weighted and unweighted) ( Figure 1B,C) and Bray-Curtis ( Figure 1D) distance, samples using a paired analysis block effects with treatment were highly variable with slight clustering by pre-and post-fertilization treatment. Adonis testing (permanova) suggested no statistically significant difference by treatment using either the UniFrac (weighted and unweighted) or Bray-Curtis difference (Table S2).  Figure S1).

Overall Metabolic Potential and Differential Metabolic Genes of Lux Arbor Switchgrass Rhizosphere
Our primary annotation enlisted KEGG for pathway and gene level metabolic potential and functions. Amongst our KEGG annotations (KO), >75% of the metabolic potential was metabolism-based (KO level 1), which was followed by Environmental Information Processing (KO level 1) at~11% (Table S3). Amongst the total metabolic potential of the Lux Arbor switchgrass rhizosphere, >15% of the metabolism-based (KO level 2) annotations were for amino acid and carbohydrate metabolism (Table S3).

Nitrogen Cycle Metabolic Potential within the Switchgrass Rhizosphere Microbiome
Nitrogen-cycling metabolic potential in switchgrass is critical to understanding why switchgrass gains little benefit from N fertilizer addition. We compared the abundance of genes related to various steps in the N cycle including N fixation (nifDHK), denitrification (nirSK, norB and nosZ), ammonification (respiration/assimilation, nrfA, napA, narG, nasA), urea catabolism (ureABC), and anammox (hzo). The urea fertilization treatment also contained urease N-(n-butyl) thiophosphoric triamide (NBPT) and nitrification dicyandiamide (DCD) inhibitors. Ammonia monooxygenase enzyme, encoded by amoA, was not found in any samples. Anammox (hzo) was not detected in any of the samples.
We further compared the N-cycling functional ORFs' impact on the whole functional metabolic potential profile using DESeq2 paired analysis, followed by MDS/PCoA ordination. Similar to the KO MDS/PCoA ordination, samples did not cluster by pre-and post-fertilization ( Figure 4A). However, as with the KO analysis, blocks 1 and 4 had greatly separated linearly post-fertilization for N-specific functional ORFs ( Figure 4A). Functionally, it appears that certain blocks were more differentiated by post-fertilization based on metabolic potential. Comparing relative abundance of N-cycling genes suggests similar abundance for ureABC, nosZ, nasA, and napA ( Figure 4B). Comparing N fixation, nifHDK gene cluster was detected in more post-fertilized samples than pre-fertilized ( Figure 4B). Blast analysis of nifD, the molybdenum-iron nitrogenase alpha chain, found that all full-length sequences were Betaproteobacterial in origin (Table S4). Of the five nifD sequences found amongst the assembled contigs, two belonged to unclassified Betaproteobacteria, two were from Dechloromonas sp., one from Sulfuriferula sp, and, lastly, one from Herbaspirillum sp. (Table S4). All nifD sequences belonged to Betaproteobacteria, with none detected from Alphaproteobacteria, even amongst the high abundance of Bradyrhizobium detected in mOTU analysis. The N fixation gene nifH, on average, had 60% higher abundance in plots pre-fertilization ( Figure 4B). Diverse members of proteobacteria contained nifHK genes, including Azonexus hydrophilus, Herbaspirillum frisingense, Herbaspirillum rubrisubalbicans, Herbaspirillum sp. HC18, Dechloromonas aromatica, Dechloromonas sp. HYN0024, Dechloromonas sp. Dech2017, and Rhodocyclaceae bacterium (Table S4). The alternative nitrogenase gene clusters anf (iron-containing nitrogenase) or vnf (vanadium-containing nitrogenase) were not detected. No rhizobium (e.g., Bradyrhizobium) nitrogenase genes were detected.
Nitric oxide reductase (norB) was detected in three out of four plots post-fertilization ( Figure 4B). NBPT had little inhibition of the gene level counts of ureABC, which encoded all the major subunits of the urease enzyme, as all samples had high levels of ureABC ( Figure 4B). Both ureABC and napA were the most abundant N cycle-related genes, regardless of the plot or timepoint ( Figure 4B). DCD also inhibits nitrous oxide (N 2 O) production when applied to the soil [61], but nosZ, which encodes nitrous oxide reductase, had similar abundance across all plots ( Figure 4B).

Genome-Resolved Metagenomics Elucidates Members of the Rare Biosphere
We obtained 190,172 > 5 kbp contigs in total from our 8 metagenomes and used the differential read abundance pre-and post-fertilization, which resulted in 29 MAGs from many phyla (e.g., Actinobacteria, Acidobacteria, Dormibacterota (formally candidate division AD3), Nitrospira, Gemmatimonadetes, Proteobacteria, and Verrucomicrobia). These MAGs represented high abundance (Actinobacteria, Acidobacteria, Proteobacteria, and Verrucomicrobia) and low abundance members (Nitrospira and Gemmatimonadetes) of common soil phyla, as well as members of the rare biosphere (Dormibacterota, Candidatus Eisenbacteria, and Candidate phyla UBA10199-formerly Deltaproteobacteria) ( Figure 6A, Table 3). Acidobacteriota (Acidobacteria) phyla had the most representative MAGs with eight ( Figure 6A, Table 3). The MAG genome sizes ranged from 2.5 to 11 Mbp, with a G + C content of 55 to 71.7% and a total contig range from 38 to 527 (Table 3). Amongst the 14 MAGs within the pre-fertilized samples, 2 are close to being high-quality drafts. Highquality drafts are defined by the minimum information about a metagenome-assembled genome (MIMAG) reporting guidelines at >90% complete, <5% contamination, and with the presence of one entire rRNA operon (5S, 16S, 23S) and 18 tRNAs [8]. The rest of the pre-fertilized MAGs were medium-quality drafts, and no low-quality drafts were used in downstream analysis. For the 15 MAGs within the post-fertilized plot, 8 were near high-quality drafts, with the rest being medium draft quality (Table 3). Table 3. Metagenomic assembled genome (MAG) assembly statistics with GTDB taxonomy. This includes the predicted taxonomy by the GTDB-Tk tool using the GTDB database, assembly statistics, completeness, contamination, and the metagenome-assembled genome (MIMAG) quality ranking. MAGs are labeled by rhizosphere soil type for pre-fertilization (I1-I14) and post-fertilization (P1-P15).    Table 3.

GTDB-Tk
We compared multiple methods to obtain taxonomic identity for MAGs, finding GTDB-Tk to be the most reliable. We classified MAG taxonomy by metawrap's classifier, classify genomes (which uses the mOTU v2 taxonomy), JSpeciesWS Tetra Correlation Search (TCS) [10,58], GTDB-Tk, and blastp of ribosomal protein S9 gene. Only GTDB-Tk provided the taxonomic identifications that were supported by ribosomal protein S9 gene blastp results. Classify genomes supplied no identifications beyond "Bacteria." (Table S5). JspeciesWS TCS misclassified candidate phyla such as Dormibacterota, which it classified incorrectly as "Mycobacterium." (Table S5). Metawraps classifier also provided no identifications for candidate phyla such as Dormibacterota (Table S5). GTDB-Tk-based taxonomy was therefore used for all downstream MAG taxonomy.
While the presence of all 29 MAGs was detected in all samples, whether pre-or post-fertilization ( Figure 6B), the abundance differed across the MAGs resolved. MAG P9 had the lowest average abundance across samples ( Figure 6B), whereas MAG P11 had the highest average overall abundance ( Figure 6B). Although we had limited replication (4 replicates), Acidobacteria and Actinobacteria phyla were also found amongst the top 10 mOTUs' phyla based on composition. The samples pre-fertilization resolved MAGs from Gemmatimonadetes, Candidatus Eisenbacteria, Nitrospira, and Dormibacterota, but these were not as highly resolved or as abundant in the post-fertilization samples (Table 3). We had one MAG in the pre-fertilization samples belonging to the "Myxococcota," formerly Deltaproteobacteria, but now its own phyla based on the GTDB. The post-fertilization samples yielded MAGs from Verrucomicrobia, Chlorflexota, and Candidate phyla UBA10199 that were not resolved well in the pre-fertilization samples (Table 3).
Gammaproteobacterial MAGs from the Lelliottia and Janthinobacterium genera (Table 3) were only found in the post-fertilization samples. MAG P6, a Lelliottia (Enterobacteriaceae), was poorly represented across all the samples except for one post-fertilization (P1) ( Figure 6B). The high abundance of P6 in a single sample could suggest an infection of plant roots within the switchgrass rhizosphere. Lelliottia are opportunistic pathogens of roots and implicated in post-harvest onion rot [63].

Betaproteobacterial MAG with Molybdenum-Based Nitrogen Fixation Gene Cluster
We screened all our MAGs for N fixation genes, such as the dominant molybdenumbased (nif gene cluster), and the alternative N fixation clusters based on vanadium (vnf gene cluster) and iron (anf gene cluster). As mentioned above, nifD genes only have similarity to Proteobacteria, with no other phyla represented. We resolved four proteobacterial MAGs in the post-fertilization treatment only. MAG P10, which is classified as Janthinobacterium, has two copies of nifD and nifK, one nifH, and one nifW, arranged on a single gene cluster for molybdenum-based N fixation. The nifHDK was previously detected as HGW-Betaproteobacteria-11 and HGW-Betaproteobacteria-7 (Table S4), which are located on the same contig within our MAG P10 genome.
Our Thermoanaerobaculia appeared to utilize nitrate but not ammonium or urea, nor to have the ability to fix N. Ammonia monooxygenase (amoA or amoB), urease (alpha or gamma, ureA or ureC), nitrification genes (nxrAB), nitrous oxide reductase (nosZ), or nitrite reductases (nirK or nirS) were not detected in the Thermoanaerobaculia MAGs. The anaerobic nitric oxide reductase transcription regulator (NorR) was present in all Thermoanaerobaculia genomes with up to nine copies in I2. Nitrate reductase (1.7.99.4; napA and nasC), formate dehydrogenase nitrate-inducible (fdnH), and nitrate transporter (narT) was present amongst the genomes. The nasA nitrate reductase was not present in any of the MAGs. MAG P4, while not classified past Acidobacteriales, had nitrate utilization genes including narT transporter, nasC nitrate reductase, and the formate dehydrogenase nitrate-inducible gene (fdnH).
We further characterized the carbohydrate utilization in the Thermoanaerobaculia MAGs for potential carbon source utilization. The most abundant genes included major glycosyltransferases (GT2 and GT4) families that synthesize diverse substrates including cellulose, chitin, sucrose, sucrose-phosphates, and glucose-glycerol phosphates. The top glycoside hydrolases (GH) encoded by the Thermoanaerobaculia MAGs included GH23 and GH0 families. GH0 is the uncharacterized family of GH, which comprises completely novel and unknown enzymes. The GH23, a specific enzyme family, is a rather specific substrate family which contains lysozyme type G (EC 3 50 (CBM50) was the most abundant CBM enzyme amongst the Thermoanaerobaculia MAGs, which contain 50 residues that have a LysM domain and works synergistically with GH23 or other enzymes that cleave chitin or peptidoglycan. CBM2 was the second-most abundant CBM enzyme present in the Thermoanaerobaculia MAGs, which has 100 residues with modules that bind cellulose, chitin, and xylan.

Nitrospira Hydrolysis of Urea, Nitrate Reduction with Limited Nitrite Reduction
Nitrospira MAGs were only well-resolved in the pre-fertilized plots. Four MAGs (I3, I8, I11, I14) in the pre-fertilized plots were classified as Nitrospira, and the genomes resolved ranged from 3 to 4.5 Mbp, 55 to 58% G + C, 82 to 96% completeness, and 4 to 8% contamination (Table 3). Nitrospira MAG I11 was the most complete at 96%, with the lowest contamination at~4% (Table 3). Nitrospira MAG I11 MAG had up to 0.11% of all the reads in a pre-fertilized plot mapped directly to the genome, representing relatively high abundance.
We further examined the N metabolism of the Nitrospira-related MAGs within Lux Arbor to identify how N metabolism functioned for these recovered MAGs. None of the MAGs had genes related to ammonia monooxygenase (amoA or amoB), so no gene annotations supported the presence of ammonia oxidation metabolism. No nitrification genes (nxrAB), nitrous oxide reductase (nosZ), or nitrite reductase (nirK or nirS) genes were detected amongst the annotations for these MAGs. All the Nitrospira MAGs had urease subunits (alpha or gamma, ureA or ureC) and accessory proteins (ureDGF) present. As for denitrification, all Nitrospira MAGs had denitrification regulatory protein (nirQ), but that was the only nir gene found amongst the MAGs. Nitric oxide reduction pathway, which encodes an anaerobic nitric oxide reductase (norV/norW), was not found; however, the anaerobic nitric oxide reductase transcription regulator (norR) was found amongst the MAGs. All MAGs had assimilatory nitrite reductase (nasE), and MAG I14 had a coppercontaining nitrite reductase, but no other genes for nitrite metabolism were found. MAG I14 had genes relating to nitrate influx and reduction to nitrite, and genes included the nitrate transporters (nasA, 3 copies) and nitrate reductase (napA, 1.7.99.4). MAG I11 had napA but did not have the nasA transporters. N fixation has never been found amongst the Nitrospira and was not found in any of our resolved MAGs.

Dormibacterota MAGs' Metabolic Potential within the Switchgrass Rhizosphere
Dormibacterota MAGs were resolved in the pre-fertilization samples only (MAGs I10 and I13). The two MAGs ranged from 3.7 to 4 Mbp, 71 to 72% G + C content, 88 to 98% completeness, and 0.9 to 2% contamination (Table 3). Dormibacterota I13 was a near highquality genome at 98.6% complete and 0.92% contamination, and was the most resolved MAG in the pre-fertilized treatments ( Table 3).
The two Dormibacterota MAGs contained no complete gene clusters for N fixation, urea, or nitrite utilization. Urea and nitrite utilization genes were also not found. I10 had a nifH nitrogenase but was missing the rest of the genes required for N fixation, such as nifDK. MAG I13 had no N fixation genes. Both had a nitrate-inducible formate dehydrogenase (fdnH) and the napA nitrate reductase. The nasC nitrate reductase and narT nitrate transporter were not found in I10 or I13 MAG. The I13 MAG lacked the nasA nitrate transporter, whereas I10 had the nasA transporter gene.
Our Dormibacterota MAGs had similar high abundance of GT2, GT4, CBM2, CBM50, GH0, GH18, GH23, and CE14 in the top 10 of their CAZy repertoire, as did our Acidobacteria and Nitrospira MAGs. This suggests that our Dormibacterota MAGs, I10 and I13, can utilize cellulose, chitin, or simple sugars such as sucrose.

Discussion
Metagenomic assembly of soil and rhizosphere ecosystems has remained an enormous challenge due to the difficulty in obtaining long contigs to reconstruct high-quality MAGs. While metagenomics has improved since the original prairie soil assembly, which contained only a few contigs > 5 kbp [9], due to the development of better software [46,47], short reads still provide significant challenges. The Lux Arbor marginal land switchgrass rhizosphere metagenomes we report represent an excellent model system with very long contigs from short reads (150 bp paired-end), due to lower microbial community diversity and complexity. In prairie soils (e.g., Kansas or Iowa) [9,10], it is not possible to obtain long contigs with just short reads and computation alone [10], and only long read technologies (i.e., moleculo) have yielded similar results to our Lux Arbor short read assembly. To compare on average, a single sample from Lux Arbor had 23,771 contigs > 5 kbp, while a Kansas native prairie soil with a similar amount of data had 4683 (100 bp paired-end) and 8532 (250 bp paired-end) [10]. A single moleculo sequence library from a pooled Kansas prairie sample yielded 10,198 contigs > 10 kbp in length, and our Lux Arbor sample yielded 44,171 contigs > 10 kbp in length using only short reads. The maximum contig length obtained from a hybrid assembly of Kansas prairie was <63 kbp, whereas Lux Arbor had 237 contigs > 100 kbp in length, with a max contig of 697,599 bp [10]. Recently, a closed bacterial genome has been obtained from the Saccharibacteria, formerly candidate phyla TM7, from stable isotope-labeled rhizosphere metagenome, suggesting that binning complete genomes directly from soil is possible [16]. Our data suggest that Lux Arbor soils have lower microbial community complexity than Kansas prairie soil, based on the quality of metagenomic de novo assembly obtained. This nominates Lux Arbor and, possibly, other marginal soils to provide a testbed for soil and rhizosphere metagenomics.
Obtaining metagenomic bins to resolve individual microbial genomes within the soil and rhizosphere has remained problematic, as the common assumption is that the higher the microbial complexity, the harder it is to resolve genomes directly from a sample. Low-complexity permafrost soil has had great success in resolving genomes, with over 1500 individual MAGs resolved with expressed metabolisms using transcriptomics and proteomics [64], but we have yet to obtain this order of magnitude with ease in nonpermafrost soil. The first genome-centric view of a soil ecosystem was in the Kansas native prairie, where 129 MAGs were obtained, but, on average, the genome completeness was quite low at~40% [10]. The second grassland soil resolved 372 total genomic bins, with 181 that were partial to near-complete [11]. A recent study of Amazon soil (using MIMAG guidelines) had 29 MAGs that were medium quality in representing over ten phyla, including members of the Candidate phyla radiation [12]. A Mediterranean grassland soil MAG study obtained 793 MAGs that were near-complete [65]. We compared concoct, maxbin2, and metabat2 within metawrap, and found that, on average, maxbin2 provided higher quality and lower contamination MAGs than metabat2, resulting in the 29 MAGs that we described ( Figure S2).
In this study, we resolved genomes in the Lux Arbor switchgrass rhizosphere that represented uncultivated phyla, including the Acidobacteria group (rare subgroup 23), Candidate phyla UBA10199, Candidatus Eisenbacteria, and Dormibacterota (AD3). Acidobacteria are dominant soil phyla, representing upwards of 20% of all soil bacteria, highly diverse, and physiologically active [66]. Acidobacteria MAGs from Kansas prairie soil were highly transcriptionally active [10], and genomes have been resolved from grassland [11] and Amazonian soil [12]. Lux Arbor Thermoanaerobaculia MAGs have previously never been described in soil, only in wastewater, sediments, and hot springs [67,68]. The Thermoanaerobaculia MAGs, we describe are the first representatives of Acidobacteria subgroup 23 from a soil or rhizosphere environment. The Acidobacteria and Dormibacterota phyla are on the 'most wanted list' of organisms from the soil and rhizosphere ecosystem for cultivation [69] and genome references via single-cell genomes or MAGs [70]. Thermoanaerobaculia and Dormibacterota MAGs have the potential to utilize nitrate, but not molecular N, urea, or ammonia. Their carbohydrate metabolism is similar to the other Lux Arbor MAGs in terms of utilization of cellulose, chitin, and simple sugars such as sucrose. The Dormibacterota have been previously implicated in carbon gas exchange (CO and CO 2 ) in Antarctic soils, including RuBisCO and carbon monoxide dehydrogenases [71]. Recently, Dormibacterota MAGs have been resolved in subsurface soil horizons which have had genes identified to aid survival in low-nutrient environments [72]. We find the Lux Arbor MAGs lack RuBisCO for CO 2 capture and utilization, but do have CO-dehydrogenases which may allow CO metabolism; CO metabolism may be thus conserved in Dormibacterota found in soil, rhizosphere, and permafrost ecosystems. Soil Dormibacterota lack both RuBisCo and dehydrogenases, meaning they lack autotrophic metabolism [72]. Dormibacterota may lose autotrophy under more stressful environmental conditions. We resolved only the third representative from the elusive candidate phyla Eisenbacteria. Our MAG is the first to be found amongst soil or rhizosphere ecosystems. The previous two were found via genome resolved metagenomics in the Atlantic Ocean deep vent sample (BioSample: SAMN09287800), named Candidatus Eisenbacteria bacterium SZUA-252, and Rifle, Colorado, USA background sediment (BioSample: SAMN04313721), named Candidatus Eisenbacteria bacterium RBG_16_71_46 [73]. This phyla appears to be extremely rare, as 8000 MAG studies did not find any representatives across thousands of samples [68]. Here, we add another representative of this rare phyla for further comparative genome analysis.
In contrast to previous studies, the metagenomes we report nominate the betaproteobacterium Janthinobacterium as a candidate organism for association nitrogen fixation in the switchgrass rhizosphere. Nitrogen fixation (nif gene cluster) was present within the bulk metagenome amongst diverse Betaproteobacterial members: Azonexus, hydrophilus, Herbaspirillum, Dechloromonas, Rhodocyclaceae, and Sulfuriferula. No other nitrogenase (nif/anf/vnf ) genes outside of betaproteobacterial class were discovered. Our reconstruction of the Janthinobacterium P10 MAG demonstrated that this genome contains a complete nif gene cluster. A related Janthinobacterium lividum V30-G6 isolated from permafrost showed low levels of N fixation via the acetylene reduction assay [74]. Bradyrhizobium spp. are highly represented in our Lux Arbor mOTU data and previous nifH data from switchgrass rhizosphere soils [40,41]. However, we were unable to find Bradyrhizobium nifDKH genes or resolve a Bradyrhizobium MAG in our study. Longer read sequencing or further depth seems to be required to address the Bradyrhizobium in Lux Arbor. Many Bradyrhizobium lack nif genes, as previously described in soils from North America and England [75,76] Two inhibitors included in the N fertilizer (SuperU) by the manufacturer inhibit nitrate reduction and urease: dicyandiamide (DCD) and N-(n-butyl) thiophosphoric triamide (NBPT), respectively. We measured the potential effects of these inhibitors at the community, individual, and gene level. Ammonium oxidation (amo) genes were not present in the bulk metagenome or associated with a resolved MAG. Dicyandiamide (DCD) limits the conversion of ammonium to nitrite via the ammonium monooxygenase (amo), but has no effect on urea hydrolysis [77,78]. DCD has been shown to sharply reduce amo gene copy numbers in ammonium oxidizing bacteria (AOB), but has relatively little effect on ammonium oxidizing archaea (AOA) amo gene copy number [78]. The AOB community is strongly shifted by DCD, with impacts on function, namely nitrification, regardless of whether it is an AOB or AOA amo [78]. This could be the reason why we were unable to detect any amo genes in post-fertilized plots, but cannot explain the lack of amo in our untreated plots. The AOB and AOA communities may be limited in abundance in these marginal lands due to lack of available ammonium that is rapidly fixed then utilized by plant roots. This could be due to rapid utilization of ammonium or its loss by volatilization, leaving little ammonium available for microbes with ammonium oxidation capabilities. In addition, lack of ammonia oxidation would lead to lower levels of nitrate available for denitrification and thus reduce potential N 2 O emissions. Indeed, DCD has been shown to mitigate N 2 O emissions [61], and other studies without DCD have observed lower N 2 O emission under switchgrass compared to other bioenergy crops. We find that the abundance of nitrous oxide reductase (nosZ) genes to be similar pre-and post-fertilization, suggesting that inhibition occurs beyond the gene level, either at an enzymatic level or due to substrate (nitrate) limitation. This is supported by previous work that found DCD did not affect nosZ gene presence or abundance [79].
Urease genes were numerically abundant both pre-and post-fertilization and were prevalent within the Nitrospira MAGs that we assembled. This points to the importance of Nitrospira in transforming urea to ammonium, and these organisms showed high abundance both pre-and post-fertilization. In another system, Nitrospira was enriched five-fold after N fertilizer treatment in agricultural soils under corn and soybean rotation, and Nitrospira MAGs were resolved with complete ammonia oxidization (comammox) [80]. Our Nitrospira MAGs lack the genes required for comammox, including the amo gene. The Nitrospira of Lux Arbor have the metabolic potential for urease, and the urease inhibitor N-(n-butyl) thiophosphoric triamide (NBPT) did not alter the ure gene copy number present pre-versus post-fertilization. However, it is unknown from our data whether urease enzyme function was altered in the post-fertilized plots, but NBPT is the most successful urease inhibitor on the market, reducing ammonium volatilization loss by 53% [81]. Further analysis of urease gene expression and urease enzyme function is needed to validate that NBPT inhibition occurs on the expression or functional level.
Carbon substrate metabolism, predicted using CAZy in the switchgrass roots of Lux Arbor, appears to be limited in terms of both the bulk metabolic potential and the individual genome level. CAZy abundances were differential under fertilizer treatment, and most were depleted in post-fertilized plots, including genes related to cellobiose, cellulose, xylan, wood-degradation, chitin, and N-acetylglucosamine. This may be due to increases in exudation of simpler carbon sources by switchgrass roots, either stimulated by fertilization or through other shifts in the system between sampling timepoints. Fertilizer treatment has been previously reported to impact CAZy enzyme function [82]. Enzyme assays could be used to validate CAZy function shifts in Lux Arbor. The MAG metabolic potentials based on CAZy had similar genes in very high abundance, with limited diversity relating to cellulose, chitin, or simple sugars.
The soils sampled here were from the same sample block, sampled just two weeks apart, immediately before and two weeks after fertilization. We found strong fertilizer effects in some blocks and little-to-no effect in others, suggesting resilience of the rhizosphere microbiome or variation in the timescales of these responses. Even with this variation, we were able to better resolve individual MAG communities present within the switchgrass rhizosphere. We are also not able to definitively conclude that the shifts we observe are due to N fertilization, since it is confounded with the time of sampling, but the functional relationships we characterize suggest that they are due to varying inputs of nitrogen and carbon in this system.

Conclusions
The Lux Arbor marginal land switchgrass plots provide an excellent model system to study a low-complexity and -diversity (<400 mOTUs) rhizosphere soil ecosystem. Due to this lower diversity, metagenomic assemblies yielded longer contigs and many MAGs at medium quality, based on MIMAG. We have described a snapshot of how a N fertilization event impacts the bulk metabolic potential of carbon and N metabolism. KO and CAZy relating to carbon and carbohydrate metabolism varied post -fertilization treatment. The two inhibitors included in the N fertilizer, DCD and NBPT, directly lowered ammonium monooxygenase (amo) gene copy number, which limits AOB but not AOA; however, our results do not confirm the impact on either functional group. AOB/AOA may have been limited in abundance in these marginal lands due to lack of available ammonium that is rapidly fixed then utilized by plant roots. We have also characterized the potential roles of several 'most wanted taxa' in the soil, resolving genomes from Candidatus Eisenbacteria, Thermoanaerobaculia, and Dormibacterota. To the best of our knowledge, our Candidatus Eisenbacteria MAG is the first described in soil at the time of this writing. Eisenbacteria are very rare in general, not being found amongst 8000 MAGs within a global survey [68]. The Thermoanaerobaculia MAGs we describe are the first representatives of Acidobacteria subgroup 23 from the rhizosphere. Our Thermoanaerobaculia MAGs utilize nitrate but not ammonium or urea, and cannot fix N. Dormibacterota have the potential for autotrophic CO utilization, which may impact carbon partitioning and storage. Furthermore, we resolved MAGs relating to ANF (Janthinobacterium) and nitrate utilization (Nitrospira). Further culture-dependent and multiomics studies are needed to evaluate the use of Janthinobacterium diazotrophs for ANF in switchgrass grown on marginal lands. Enhancing microbial ANF within marginal land sites will decrease fertilizer usage, lowering N needs for biofuel productions. Using ANF over industrial N fertilizer will require less energy overall, produce less carbon dioxide emissions, and provide a potential greener future of energy.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/agronomy13051294/s1, Figure S1: Field site map with sampling blocks in Lux Arbor. Figure S2: Comparison of metabat2 and maxbin2 binning statistics. Table S1: mOTU adonis statistical testing within phyloseq R results. Paired analysis using sample block (B1-B4) "Block," and post-fertilization. Table S2: Top-10 mOTU abundance table from phyloseq R. This includes taxonomy of mOTUs and sample metadata. Table S3: KEGG KO contig proteincoding ORF count table. This includes level I to level II raw counts for KEGG KO. The DESeq2 KO paired analysis significant table out of 3204 KO's with nonzero total orf count (p-value < 0.05). Table S4: Nitrogenase molybdenum-iron protein alpha chain gene (nifD) blast table. This was blastp analysis against genbank/refseq, with score, accession, and taxonomy. Table S5: Comparison of MAG taxonomy annotation. Comparing metawraps MAG taxonomy tool, Classify genomes tool, Jspecies, and ribosomal gene S9 with genbank/refseq, with score, accession, and taxonomy. Data Availability Statement: Raw sequence data, assembled contigs, and supplemental data are all available on https://osf.io/mzrvj/ (accessed on 15 April 2023). All code for this study is available on www.github.com/friesenlab/MMPRNT_panicum_metagenome_mags/ (accessed on 15 April 2023).