Aboveground and belowground responses to cyanobacterial biofertilizer supplement in a semi‐arid, perennial bioenergy cropping system

The need for sustainable agricultural practices to meet the food, feed, and fuel demands of a growing global population while reducing detrimental environmental impacts has driven research in multi‐faceted approaches to agricultural sustainability. Perennial cropping systems and microbial biofertilizer supplements are two emerging strategies to increase agricultural sustainability that are studied in tandem for the first time in this study. During the establishment phase of a perennial switchgrass stand in SW Montana, USA, we supplemented synthetic fertilization with a nitrogen‐fixing cyanobacterial biofertilizer (CBF) and were able to maintain aboveground crop productivity in comparison to a synthetic only (urea) fertilizer treatment. Soil chemical analysis conducted at the end of the growing season revealed that late‐season nitrogen availability in CBF‐supplemented field plots increased relative to urea‐only plots. High‐throughput sequencing of bacterial/archaeal and fungal communities suggested fine‐scale responses of the microbial community and sensitivity to fertilization among arbuscular mycorrhizal fungi, Planctomycetes, Proteobacteria, and Actinobacteria. Given their critical role in plant productivity and soil nutrient cycling, soil microbiome monitoring is vital to understand the impacts of implementation of alternative agricultural practices on soil health.


| INTRODUCTION
Over the last century, conventional agricultural practices, including synthetic fertilization of nitrogen (N), tillage, and the use of chemical herbicides, have benefited humanity by drastically improving agricultural productivity, but at a significant cost to soil and ecosystem health (Chen et al., 2014;Rockström et al., 2009). For example, agricultural intensification and mismanagement of fertilizer application have led to soil acidification, decreased soil organic carbon (SOC), nutrient loss through leaching and gaseous emission, and declines in microbial diversity and functioning (Barak et al., 1997;Tian et al., 2020;Tsiafouli et al., 2015). The degradation of soil microbial communities in particular has become an increasing concern as microbes including bacteria, archaea, and fungi are critical regulators of soil organic matter dynamics, nutrient availability, and subsequent crop productivity (Saleem et al., 2019;Singh et al., 2020;van der Heijden et al., 2008;Wagg et al., 2014). To mitigate or reverse these effects, the expansion in both research and practice of sustainable cropping systems has become a focus to address negative environmental impacts while meeting, or surpassing, conventional cropping system yields (Hobbs et al., 2008). The use of perennial cropping systems and the use of microbial biofertilizers as supplements to synthetic fertilizers are two emerging strategies to address the negative impacts of agriculture on the environment, and specifically the mismanagement of nitrogen (N) fertilizer (Robertson & Vitousek, 2009). These approaches provide multiple key ecosystem services (improved C cycling and microbial diversity, reduced nutrient loss) in comparison to conventional agricultural practices (Jesus et al., 2016;Robertson & Vitousek, 2009;Werling et al., 2014), and in combination could result in dramatically improved sustainability.
Perennial crop cultivation has increased significantly in the United States in recent years, particularly in the use of switchgrass (Panicum virgatum) for cellulosic bioenergy production Robertson et al., 2017;Wright & Turhollow, 2010). Switchgrass cropping systems may also play an integral role with emerging climate-smart technologies such as bioenergy with carbon capture and storage (BECCS) because deep-rooting perennial plants transfer large amounts of carbon belowground (Jarchow et al., 2015;Stoy et al., 2018). While it remains unlikely that any single species will become a universal feedstock for bioenergy production (Jessup, 2009), perennial grasses have many environmental advantages (enhanced biodiversity, pest suppression, and pollination) compared to annual species (de Oliveira et al., 2018;Fazio & Monti, 2011;Lai et al., 2018;Pugesgaard et al., 2015;Werling et al., 2014). Furthermore, there is high economic potential for some perennial bioenergy crops (Adler et al., 2007;Wright & Turhollow, 2010), and although spatial and temporal variations in yield can be high (Schmer et al., 2009), consideration of regional biophysical factors during species selection can help optimize production. However, the net benefits of conversion to perennial bioenergy systems depend largely on the retention of ecosystem carbon (C) and N (Robertson et al., 2011;Robertson & Vitousek, 2009). While perennial grass systems tend to accumulate soil organic carbon (SOC) when converted from annual-based systems (Lemus & Lal, 2005;Liebig et al., 2008), variability in SOC accrual has been documented across agroecosystems (Ye & Hall, 2019). In modeling scenarios developed by Dolan et al. (2020), future climate projections indicated that much of the Upper Missouri River Basin (UMRB) of North America has high potential for second-generation bioenergy production. However, the productivity of switchgrass, a C4 species adapted to warm climate regions (Sanderson et al., 2006), has not been rigorously documented in field experiments across agricultural regions of the semi-arid UMRB. Given the extent of marginal land and relatively low quality of soils across the UMRB, the ecosystem services provided by perennial cropping systems could be highly beneficial (Emery et al., 2017;Gelfand et al., 2013). As in all bioenergy cropping systems, nutrient, and specifically N, loss in annually harvested biomass must be matched with appropriate annual inputs. Thus, combining a perennial cropping system with a biofertilizer might offer a sustainable approach to replacing annual N loss (from harvest) through the input of an initially less reactive form of N that is not as likely to be lost quickly from volatilization or leaching.
Biofertilizers are inoculants of living microorganisms applied to plants or soil to support plant growth by increasing nutrient availability (Bhattacharjee & Dey, 2014;Pellegrino et al., 2015). N-fixing cyanobacteria are gram-negative bacteria that often exhibit filamentous growth and are widespread in soil and aquatic habitats, and that have been used as biofertilizers for many years in rice paddy operations (Alvarez et al., 2021;De, 1939;Prasanna et al., 2015). They also can promote the growth of a wide variety of terrestrial crops with additional benefits to soil quality (Nain et al., 2010;Swarnalakshmi et al., 2013;Watanabe et al., 1979). Cyanobacterial-based biofertilizers (CBF) have been shown to (1) maintain plant productivity when used as a supplement to replace 25%-50% synthetic NPK fertilizer (Chittapun et al., 2018;El-Beltagy et al., 2016), (2) increase SOC through photosynthesis and biomass turnover (Yilmaz & Sonmez, 2017), (3) decrease soil erosion through the production of exopolysaccharides (Falchini et al., 1996;Malam Issa et al., 2006), and (4) increase a variety of soil enzymatic activities (Nisha et al., 2018). In addition, when compared to synthetic fertilizers, CBF provides a more gradual release of nutrients to the soil which can benefit crop growth and quality (Coppens et al., 2016;Mulbry et al., 2007). In contrast, there is ample evidence that urea application, which accounts for >55% of global N production and is the most dominant crop N fertilizer (FAO, 2019), is susceptible to large losses of ammonia (Pan et al., 2016) and nitrous oxide (Bouwman et al., 2002) to the atmosphere and nitrate to groundwater (Di & Cameron, 2002). In regards to soil N dynamics in bioenergy cropping systems, if fertilization of N is managed to match crop needs (Ruan et al., 2016), perennial systems have a greater potential for N retention in comparison to annual crops via reduced loss of nutrients through leaching, runoff, and greenhouse gas flux coupled with increased N use efficiency (Culman et al., 2013;Liebig et al., 2008). While the reported benefits of CBF are numerous, the effects of CBF on the soil microbial community and nutrient cycling in field studies remain a large knowledge gap. The vital roles soil microbes play in decomposition and nutrient cycling make understanding the impact of biofertilizers on native soil microbial communities critical to their successful implementation in sustainable agricultural systems.
Together, perennial bioenergy crops and biofertilization could significantly promote agricultural sustainability by improving nutrient retention and microbially mediated ecosystem processes. To our knowledge, this is the first study to investigate these strategies in tandem. The goal of this study was to analyze aboveground and belowground treatment effects of CBF on perennial crop establishment, soil microbial community diversity, N availability, and SOC in a field-based, randomized plot experiment. We hypothesized that (1) supplementing synthetic fertilizer with CBF would maintain aboveground biomass production compared to synthetic fertilization alone and that the slower mineralization rate of CBF would benefit lateseason crop productivity, (2) CBF supplementation would prolong N availability, and (3) CBF supplementation would benefit microbial diversity compared to synthetic fertilization alone and response to fertilization would be phylogenetically conserved in bacterial/archaeal and fungal communities.
To test these hypotheses, we compared aboveground and belowground responses of switchgrass field plots treated with synthetic urea fertilizer, a mixture of urea and CBF, or untreated (controls). As perennial grasses often take a year to establish and can suffer from weed competition in the early establishment phases, barley was intercropped to outcompete weeds and served as an early indicator of fertilizer effects, as barley has a higher N demand than switchgrass (Ball & O'Sullivan, 1987;Lemus et al., 2008). To preserve the switchgrass stand integrity during establishment, we measured aboveground productivity non-destructively using the Normalized Difference Vegetation Index (NDVI) as a proxy for aboveground biomass (Garroutte et al., 2016). Soil chemistry and microbial diversity using high-throughput sequencing were measured throughout the growing season to determine the chemical and biological effects of CBF implementation.

| Site description and experimental design
Field plots are located at the Montana State University (MSU) Arthur H. Post Research Farm (45.66 N, −111.15 W) located 10 km west of Bozeman at 1450 m elevation where mean annual temperature and precipitation are 6.5°C and 408 mm year −1 ( Figure S1). The soil at the Post Farm is broadly classified as an Amsterdam silt loam which is a fine-silty, mixed, superactive, frigid Typic Haplustoll (100 g kg −1 clay, 810 g kg −1 silt, 90 g kg −1 sand; Engel et al., 2017). The preceding crop was Proso millet (Panicum milaceum) which was fertilized with 67 kg N ha −1 urea in spring 2017. Barley (Hordeum vulgare, Hayse forage variety) and switchgrass (Panicum virgatum, Dakota variety) were seeded on April 16, 2018 at 16 kg ha −1 with a no-till disk seeder. Three treatments were established across replicated (n = 5) 3.34 m 2 plots: (1) 100 kg N ha −1 urea (urea), (2) 100 kg N ha −1 50:50 CBF/ urea mix (CBF:Urea), and (3) control plots that received no experimental N addition (control). To minimize edge effects, the plots were arranged in a randomized block design with buffer alleyways (1.2 m) between replication rows that were seeded but left untreated ( Figure S2). We elected to supplement the synthetic fertilizer with a halfdose of CBF as recent literature found replacing 50% of the synthetic fertilizer application with CBF resulted in the greatest yields of various crops compared to synthetic fertilizer alone (Chittapun et al., 2018;Naher, 2018;Nain et al., 2010). Notably, we were limited by production capacity to expand the field plot treatments to include a CBF-only treatment to separate CBF effects from urea effects. However, this design allowed us to build off of the existing literature and study CBF in a practical manner similar to what producers would typically utilize (Naher, 2018). Fertilizer application timing was based off common regional application practice for urea and was applied on April 22nd 2018 (Jacobsen et al., 2005). Urea plots were hand-spread and received all 100 kg ha −1 in April, whereas CBF:Urea plots only received 50 kg ha −1 . The remaining 50 kg ha −1 of CBF were split into five application periods (June 2, June 15, June 28, July 12, and August 09) applied in 10 kg ha −1 increments explained in detail below. Supplemental sprinkler irrigation (n = 4 events), totaling ~112 mm of precipitation equivalent, was applied after fertilization applications to all plots. This was done to help the biofertilizer infiltrate below the soil surface when natural precipitation events were lacking. Due to cost constraints, soil chemistry and microbial community results are reported from a subset of 3 out of the 5 treatment replicates.

| CBF cultivation
From previous screening of several cultured representatives of N-fixing cyanobacterial species, we selected Anabaena cylindrica (UTEX 1611) for use as CBF as we verified that it grows robustly in liquid media, has a relatively high N fixation rate in liquid culture, and does not produce any known cyanotoxins (Weeks, 2013). Elemental analysis was conducted on an EA Flash 2000 and indicated that CBF biomass contained 9% N, 1% P, 1% K, and 50% C. For each application period, we cultivated a total of 800 L liquid culture using 3 × 200 L raceway ponds and 4 × 50 L bag reactors (Figure 1). At the time of harvest, bags and raceways were pumped into conical tanks to gravity settle overnight before harvesting the concentrated biomass. Dry biomass weights (in mg/L) were measured in triplicate by vacuuming the concentrated culture onto 0.2 μM filters and drying overnight at 60°C. Dry weights were used to calculate the volume of culture to be applied to each plot. CBF was cultured for 2 weeks prior to harvesting, resulting in five CBF applications over the growing season. Due to limited production capacity, living CBF biomass was supplemented with commercially available, dried Spirulina (Arthrospira platensis) biomass (10% N, 1% P, 1% K) at a final biomass ratio of 20:80 A. cylindrica:Spirulina to meet the goal of 50% CBF-derived N for the CBF:Urea supplement treatment. The multiple smaller applications and mixture with Spirulina were necessary due to limitations in culture volume of A. cylindrica.

| Sampling and analyses
Crop productivity and growth were monitored in situ using a novel remote sensing device (Arable, USA) designed to measure a suite of crop management indices including hourly spectrometer data. A wide variety of vegetation indices are used to assess vegetation/crop status via remote sensing devices metrics. Here, we use the Normalized Difference Vegetation Index (NDVI), one of the most widely used vegetation index, that has been correlated to aboveground biomass, canopy area, and leaf N content (Cabrera-Bosquet et al., 2011;Teal et al., 2006;Tucker, 1979), where NDVI = (NIR − red)/(NIR + red) and NIR is the fraction of near-infrared radiation that is emitted back to the sensor and red is the fraction of red radiation emitted back to the sensor. The index reports values from −1 to 1 and generally represent the vegetative greenness within the spectrometer's footprint. Arable Mark devices were installed directly following seeding on April 26, 2018. Arable Marks were mounted on 2 m tall steel posts (7.5 cm diameter) in the center of each plot, installed at 1 m above maximum crop height, and stayed in place for the duration of the growing season ( Figure  S3). In all, six Arable Marks were deployed for the study across two randomly selected plots from each treatment. The raw mean daily NDVI values from each Arable sensor were aggregated by month to test for statistical differences between treatments across the growing season. Barley biomass was harvested July 23, 2018 at peak forage biomass by shearing 10 cm above ground surface in 1 m 2 blocks at the center of each plot. The biomass was air-dried for 1 week at 40°C and subsequently weighed as dry biomass.
Field plots were sampled prior to seeding in April 2018 to establish background soil conditions and then monthly to capture temporal and compositional variance in the soil microbial community and shifts in N availability. Due to the small area of individual plots, one fixed-depth soil core (Oakfield Apparatus) was taken from each plot and broken up into 0-15 and 15-30 cm soil depth increments. Three plot replicates were haphazardly selected for soil chemical and microbial community analyses. The soil was then homogenized in separate pools by depth and each pool was subsampled for microbial and chemical analyses. The biological subsamples were immediately flash frozen in an ethanol-dry ice bath and temporarily stored on dry ice before transfer to long-term storage at −80°C. Before soil processing, a small subsample was separated and massed before being dried at 105°C for 48 h and massed again for volumetric soil moisture content ( Figure S4). To determine bulk soil total nitrogen (TN) and soil organic carbon (SOC), samples were immediately transported back to the lab and oven-dried at 65°C for at least 24 h, passed through a 2 mm sieve, and milled for 12 h in stainless steel jars on a F I G U R E 1 Cyanobacterial culturing systems in the Montana State University Plant Growth Center greenhouse: 200 L raceway ponds (left) and 50 L photobioreactors (right) slow rolling tumbler. Samples were prepared in duplicate as SOC quantification required acidification of the soil using 1 and 2 M HCl to remove the inorganic C fraction. After drying for 48 h, %C and %N were determined using a Costech ESC 4010 elemental analyzer. Duplicates were run for ~10% of all samples to ensure analytical precision. Inorganic nitrogen (NO 3 − -N and NH 4 + -N) was processed within 12 h of soil collection by using ~10 g of the homogenized field moist soil core, with coarse/fine root material removed, and then extracted in 100 ml of 2 M KCl by shaking at 200 RPM for 30 min. Following shaking, the samples were filtered through Whatman grade-2 filter paper and frozen until analysis. We determined extractant concentrations of NO 3 − -N and NH 4 + -N colorimetrically using an automated flow injection system (Lachat QuickChem 8500, Hach Company). NH 4 + -N was determined using the indophenol method while NO 3 − -N was determined using a cadmium column reduction and coupled colorimetry. Soil pH was measured in 2:1 deionized water to soil solution using a glass electrode (Mettler Toledo).

| DNA extraction, amplification, and sequencing
To obtain high-resolution representation of the soil microbial community and its response to fertilization, we performed rRNA gene amplicon sequencing across three of the treatment replicates at three time points (n = 54). The selected cores were specifically paired with those used for inorganic-N extraction to directly link measures of the microbiome and N availability. Genomic DNA was extracted from 0.5 g frozen soil samples using the Fast DNA SPIN Kit for soil (MP Biomedicals) according to the manufacturer's instructions. We targeted the V4 region of the 16S rRNA gene using f515a and r806b (GTGYCAGCMGCCGCGGTAA and GGACTACVSGGGTATCTAAT, respectively) from the Earth Microbiome Project (Thompson et al., 2017), and the D2 hypervariable region of the fungal LSU using LR22r and LR3 (Mueller et al., 2016). PCRs were carried out in 20 μl volumes with each standard reaction mix containing 1 U Phusion high-fidelity polymerase (ThermoFisher Scientific), a final concentration of 1x Phusion HF reaction buffer, 200 μM dNTPs, 0.5 μM each primer, 1 μM BSA, 1 mM MgCl, sterile molecular grade water, and 5 μl (10 ng) genomic DNA as a template. PCR conditions consisted of initial denaturation at 98°C for 30 s, followed by 22 cycles of 98°C for 30 s, 55°C for 30 s, 72°C for 45 s and a final extension at 72°C for 5 min. Barcoding was conducted using the Illumina Nextera Indexing Kit D, with eight cycles of indexing PCR in 20 μl volumes using the same concentrations as PCR1 with 10 μl template DNA.
PCR products were checked for quality and length in a 1% agarose gel and quantitated using the Quant-iT™ dsDNA Kit (Invitrogen) with a BioTek H2 plate reader. The two libraries were pooled at 4 nM and loaded onto an in-house Illumina MiSeq (Illumina). Reads were merged, trimmed, and dereplicated with USEARCH, and zero-radius operational taxonomic units (ZOTUs) were identified with UNOISE3 (v.11.0.667;Edgar, 2010). Similar to amplicon sequence variants, UNOISE3 identifies ZOTUs by correcting point errors and filtering out chimeric amplicons to obtain accurate predictions of true biological sequences that are superior to utilizing a 97% identity cutoff (Callahan et al., 2017). Representative 16S and LSU ZOTUs were classified against their respective databases using the online Bayesian classifier of the Ribosomal Database Project (Set 18 for bacteria, Set 11 for fungi; Wang et al., 2007). All reads classified to chloroplast (16S), or animalia, protozoa, or viridiplantae (LSU) were removed prior to analysis, and datasets were randomly subsampled according to the sample with the lowest read number for statistical analysis using the R package vegan (Oksanen et al., 2019). The data are accessible at MG-RAST under accession numbers: mgm4916997.3 (LSU) and mgm4916998.3 (16S).

| Phylogenetic trees
Reference phylogenies for the 16S and fungal LSU were constructed using full-length and near full-length sequences downloaded from GenBank. For fungi, we focused on sequences associated with the Assembling the Fungal Tree of Life project (Kauff et al., 2007). For bacteria and archaea, we used sequences from both isolates and metagenome assemblies to better cover the tree of life, including a subset of the Candidate Phylum Radiation. Reference sequences were aligned using mafft (Katoh & Standley, 2013;Stamatakis, 2006) and maximum likelihood trees were constructed using RAxML with the GTR + gamma model. ZOTUs were aligned and placed onto their respective reference phylogenies using pplacer (Matsen et al., 2010). Trees were visualized and annotated using the interactive tree of life (iTOL; Letunic & Bork, 2019).

| Data analysis
To test the effect of fertilizer treatment on soil chemistry by month, we used linear mixed models using fertilizer treatment and depth as fixed effects and plot replicates as a random effect (to account for the non-independence of replicate soil samples taken within a replicate) using the lme4 package (Bates et al., 2015). To assess whether there were differences in soil chemistry based on fertilizer treatment, we performed type III ANOVA; if fixed effects accounted for variation, we compared post-hoc differences using pairwise Tukey comparisons using the emmeans package (Lenth et al., 2020). For NDVI analysis, we fit linear models with a fixed effect of treatment and fit one-way ANOVA's. If treatment accounted for variation, we used Tukey's comparison to examine pairwise comparisons for fertilizer treatment effects. We report all soil chemistry and NDVI results as model mean values (±1 SE) and C:N as molar ratios. For soil dissolved inorganic nitrogen (DIN; NO 3 − -N and NH 4 + -N) analysis, we compared concentrations at the individual soil depths (0-15 and 15-30 cm) as well as the mean aggregated total depth (0-30 cm). pH is reported as the mean of hydrogen ion concentrations followed by a log10 back transformation with uncertainty reported as a 95% confidence interval (CI). For NO 3 − -N and NH 4 + -N samples that were under the detection limit, we substituted a value that was ½ the lower detection limit. All data were assessed for influential points, homogeneity of variance, and normality. Log transformations were performed on data that did not fit the assumption of normality. Treatment differences were declared significant at p = .05 for NDVI and microbial community analyses and p = .10 for NO 3 − and NH 4 + because of the inherent temporal and spatial variability of soil N transformations and the limited number of cores taken per plot (Moebius-Clune et al., 2008).
We calculated alpha diversity and beta diversity of the bacterial/archaeal and fungal communities using the phyloseq package (McMurdie & Holmes, 2013) with three separate methods for alpha diversity: ZOTU richness, Shannon's diversity index (Shannon, 1949), and Faith's Phylogenetic Diversity (PD;Faith, 1993), and unweighted UniFrac for phylogenetic-based beta diversity (Lozupone & Knight, 2005). We used PERMANOVA to test for treatment, soil depth, and seasonal effects on community beta diversity (Anderson, 2017). To assess individual ZOTU responses to fertilizer treatment, response ratios were calculated by taking the difference between the read count for each ZOTU in a fertilized treatment and that of the control, non-fertilized field plots divided by the total number of reads for that ZOTU in both the treatment and control: The response ratio calculation was permuted 999 times for each ZOTU and averaged to acquire a final response ratio. The output of Equation (1) is the responses on the scale from −1 to 1. Response classifications were determined by splitting the responses into quartiles. To investigate differences in response ratios on the scale of individual ZOTUs, the vegan package was utilized to conduct SIMPER analysis of ZOTUs driving the differences between the response ratios in treatment groups (Oksanen et al., 2019). The phylosig package was utilized to test Pagel's lambda and Abouheif's C mean measures of phylogenetic signal (Keck et al., 2016). Prior to calculation, ZOTU tables were trimmed to include only ZOTUs with 10 or greater total reads to reduce erroneous classifications due to sampling error. We performed all analyses in R version 3.6.1 (R Core Team, 2019).

| Barley and switchgrass growth
Fertilizer treatment influenced NDVI throughout the growing season (F 10,952 = 5.37, p < .0001). Compared to the control and urea treatments, mean daily NDVI values were higher in the CBF:Urea treatments during the second half of the 2018 growing season (Figure 2; Table  S1). Barley plant biomass was similar between the urea and CBF:Urea treatments (p = .189) with 8806 ± 1362 and 7288 ± 774 kg ha −1 , respectively ( Figure S5). However, both urea and CBF:Urea treatments significantly increased aboveground plant biomass compared to the unfertilized control at 4170 ± 1143 kg ha −1 (p = .001 and p = .02, respectively).
F I G U R E 2 Box plots (median, 5th and 95th percentiles) of daily mean NDVI values (n = 980) aggregated by month across the three fertilizer treatments. The dashed line represents when the barley intercrop was harvested. Differences (p < .05) between aggregated daily means (Table S1) are denoted by letter signifiers ratios at the 0-15 cm soil depth declined by 16% and 10% across all treatments (F 1,14 = 22.88, p = .0029, F 1,14 = 6.02, p = .025) from pretreatment (April) to the end of the season (September). Total soil N and soil pH did not vary across fertilizer treatments or through time (Tables S2 and S3). Experimental treatments had two major effects on soil DIN across the cumulative growing season (May-September). First, the CBF:Urea treatment resulted in overall higher mean concentrations of NH 4 + -N compared to unfertilized controls (0-30 cm, p = .01), whereas the urea treatment concentrations were not different than the unfertilized controls (0-30 cm, p = .36). There was suggestive evidence that this effect was stronger at the 0-15 cm depth (p = .09). Second, the urea treatments had higher mean concentrations of NO 3 − -N at the 15-30 cm depth compared to the unfertilized controls (p = .05) while there was no difference between the CBF:Urea and control treatments (p = .43).
Additionally, we found evidence for seasonal influence over fertilizer treatment effect on DIN (F 4,139 = 13.95, p ≤ .0001; Figure 3). During the spring (May and June), fertilizer treatment affected the NO 3 − -N concentrations at both soil depths (F 2,54 = 12.3, p ≤ .001, F 2,54 = 3.37, p = .071; Figure 3). Urea concentrations were greater than the unfertilized control (p = .006, Figure 3) at the 15-30 cm depth while there was no difference between the CBF:Urea and control treatments (p = .16, Figure  3). In contrast, when just analyzing the spring months (May and June) NH 4 + -N concentrations, we found no differences among fertilizer treatment or depths. During the late season (August and September), we report evidence that fertilizer treatment affected DIN concentrations (F 2,53 = 5.11, p = .009). The CBF:Urea treatment had higher concentrations of NH 4 + -N compared to the unfertilized controls (p = .007; Figure 3) with marginal evidence of a difference between CBF:Urea and urea treatments (p = .12). Furthermore, CBF:Urea was the only treatment with a difference in NH 4 + -N concentrations between the two soil depths (p = .02), with higher concentrations at 0-15 cm. In contrast, when analyzing just the late season, we did not find any differences in NO 3 − -N concentrations among fertilizer treatment or soil depths.  Figure 4a). Composition of the bacterial/archaeal community at the phylum level was stable across time, depths, and treatments ( Figure 4a). The fungal community comprised Ascomycota (38%), Basidiomycota (21%), Glomeromycota (21%), Chytridiomycota (9.7%), Blastocladiomycota (2.9%), and the remaining were unclassified fungi incertae sedis (6.21%; Figure 4b). The Glomeromycota was the most dynamic phylum which increased in relative abundance over time from 13% in April to 32% in September while the Ascomycota decreased from 44% to 30% in this timeframe. In addition, Glomeromycota had a greater relative abundance in the deeper soil in general (30% vs. 15%) and was also slightly greater in the urea-treated plots than in the CBF:Urea and control plots (25% vs. 18% and 22%, respectively).

| Cyanobacterial response to fertilizer treatment
Cyanobacteria were rare, averaging only 0.07% of the total community over the entire season. At the genus level, the cyanobacteria were dominated by Microcoleus sp. and Nodosilinea sp., and BLASTn of these ZOTUs indicated that Microcoleus vaginatus was the primary species present ( Figure 5). The CBF species Anabaena cylindrica was not detected in any field plots prior to application and was only detected in CBF:Urea plots in July.

| Bacterial/archaeal response to fertilizer treatment
As shown in Figure 6, alpha diversity of the bacterial/ archaeal community had similar trends among diversity metrics with no strong evidence for differences by treatment, time, or soil depth except for phylogenetic diversity (Figure 6a,b). There was strong evidence of a seasonal shift in phylogenetic alpha diversity with an increase from April to July with community mean branch lengths of 802 ± 18 and 875 ± 16 (SE), respectively (F (2,48) = 5.41, p = .008; Figure 6b). September PD was similar July PD with a mean branch length of 853 ± 14. CBF:Urea plots trended toward greater alpha diversity although there was high variability between plots. The beta diversity of the bacterial/archaeal community indicated a strong separation by soil depth (F (2,48) = 3.49, p = .003) but not by treatment (p = .939) or time (p = .088; Figure S6a, Table  S6), suggesting that distinct communities were present at each soil depth, but the respective levels of diversity were similar.
Response ratios of the changes in relative abundance of ZOTUs across the phylogenetic tree indicated that bacteria generally had neutral responses across the phylogenetic tree ( Figure 7a). However, Pagel's λ indicated a significant phylogenetic signal in the bacterial/archaeal (λ = 0.23, p = 2.60e −101 ) and confirmed with Abouheif's C mean (Table  S7) meaning that response ratios of closely related ZOTUs were more similar to each other than would be expected under random trait distribution. Density plots of response ratios confirmed that the bacterial/archaeal community had primarily neutral to weak responses (Figure 7a insert) and that the differences in responses between the two fertilizer treatments were slight. However, the CBF:Urea treatment had slightly more strong negative and weakto-moderately strong positive responses than the urea treatment while the urea treatment had more neutral responses. Through our SIMPER analysis, we identified 12 ZOTUs driving the differences in response ratios between the urea and CBF:Urea bacterial/archaeal community (Table S8; p < .001). The four most influential ZOTUs contributing to were classified to Planctomycetaceae and two Actinobacterial families, Micromonosporaceae and Acidimicobiaceae while the others also included Proteobacteria, Candidatus Uhrbacteria, Acidobacteria, Tenericutes, Abditibacteriota, and Bacteroidetes.

F I G U R E 4
Relative abundance at phylum level of (a) bacterial/archaeal and (b) fungal community phyla by treatment, soil depth, and time F I G U R E 5 Relative abundance of cyanobacterial genera across time, treatments, and soil depth 3.6 | Fungal response to fertilizer treatment Similar to the findings for the archaea/bacteria, we also observed also high variability in the alpha diversity of the fungal community with weak evidence for a difference between treatments regardless of the diversity metric. However, there was evidence for a difference in taxonomic richness over time (F (2,48) = 5.59, p = .007; Figure  6c). The number of fungal ZOTUs was highest in April (274 ± 11), followed by September (250 ± 9) and July (225 ± 13). Like the bacteria, beta diversity of the fungal community indicated a strong separation by depth but not by treatment or time ( Figure S6b). Through SIMPER analysis, we identified 15 ZOTUs driving the difference between the urea and CBF:Urea fungal communities, which appeared to be primarily driven by arbuscular mycorrhizal fungi (AMF). The four most influential ZOTUs driving 76% of the difference between the two communities all belonged to the Glomeromycota, with the urea treatment having greater abundance than the CBF:Urea treatment (Table S9). The remaining ZOTUs were classified as Eurotiomycetes, Monoblepharidomycetes, Agaricomycetes, Pezizomycetes, Sordariomycetes, and Zygomycota incertae sedis. However, the two communities were found to be 74% different overall, with individual ZOTUs making up at most 0.07% of this difference.
In contrast to the bacterial/archaeal community, fungi exhibited polarized responses to fertilization with primarily strong positive or negative responses from either urea or CBF:Urea treatments compared to the non-fertilized control (Figure 7b). Phylogenetic signal was also detected in the fungal community (λ = 0.621, p = 9.9e −37 ), which was confirmed with Abouheif's C mean (Table S7) suggesting that, albeit weakly, response ratios of closely related ZOTUs were more similar to each other than would be expected under random trait distribution. The fungal community was polarized with primarily strong positive or strong negative responses and confirmed the similarity between the two treatments (Figure 7b insert). In contrast to the bacterial/archaeal community, the urea treatment had slightly more strong negative and strong positive responses in the fungal community than the CBF:Urea treatment and the CBF:Urea treatment had a greater frequency of moderate-to-neutral responses.

| DISCUSSION
Maintaining crop productivity while reducing negative environmental impacts is critical to increasing the sustainability of agricultural practices. Perennial cropping systems may play an important role in future sustainability initiatives, but their implementation in dryland agriculture is limited and as a result the impacts on soil C and N dynamics are unclear (Wang et al., 2015). In addition, understanding the complexities of how the soil microbiome effects ecosystem function and contributes to soil health is currently a top research priority in the field of microbial ecology (Bender et al., 2016;Fierer et al., 2007;Mendes et al., 2013). Furthermore, field experiments and data demonstrating the effects of algal biofertilizers on the structure and functioning of microbial communities in agricultural systems are scarce. Therefore, this study aimed to contribute to these knowledge gaps by evaluating the effects of CBF supplementation on switchgrass stand establishment and to characterize the initial soil microbial community and N availability.

| Fertilizer effects on plant productivity
Supplementing 50% of the synthetic fertilizer application with CBF compared to 100% urea application resulted in similar barley biomass and provided field evidence for sustained increases in late season NDVI for an establishing perennial switchgrass crop (Figure 2). Similarly, in a study with potted rice, a combination of Nostoc sp. CBF with a half-dose of synthetic fertilizer resulted in the highest yield of rice grain compared to CBF or synthetic fertilizer alone (Chittapun et al., 2018). It is likely that by combining synthetic and biofertilizers, plants can take advantage of readily available synthetic N then utilize the F I G U R E 6 Alpha diversity of bacterial/archeal (a, b) and fungal (c, d) communities by treatment and time. Two measures of alpha diversity, #OTUs (A, C), and phylogenetic diversity (mean branch length), (b, d). Error bars represent standard error (SE) from the mean (n = 6). Statistics compare diversity between time points based on Tukey-Kramer multiple comparison tests. *p < .05 and **p < .01; ns, non-significant more slowly mineralized biofertilizer biomass throughout the growing season (Coppens et al., 2016;Kubheka et al., 2020).

| N and C responses to fertilization
Studies carried out with various CBF species have often shown that CBF can maintain or improve upon crop growth compared to synthetic fertilizers with additional benefits to soil quality including increased nutrient availability longer than synthetic fertilizers (Manjunath et al., 2016;Prasanna et al., 2012). We observed significantly higher late-season NH 4 + -N concentrations in shallow soil across the CBF:Urea plots compared to the urea fertilized plots or the non-fertilized controls in August (Figure 3). Furthermore, Anabaena CBF was detectable in July but not in September ( Figure 5), which when considered with the increased soil NH 4 + in the CBF:Urea plots, suggesting that the algal biomass was mineralized into available N over that time period, as observed by Coppens et al. (2016). This highlights the potential for biofertilizer supplements to extend N availability well into the late growing season and support a more closed N cycle in comparison to the conventional N application via urea pellets. We also observed higher levels of NO 3 − -N in the urea-only treatment at the 15-30 cm soil depth across the entire growing season, and especially in the early season ( Figure 3). This is noteworthy as NO 3 − -N is a highly soluble and mobile form of N that, at this soil depth and early stage in stand development, is highly susceptible to eventual groundwater loss. We did not expect to see shifts in bulk soil chemistry after one growing season (Schrumpf et al., 2011); however, the ubiquitous decline in soil SOC across all treatments (Table S2) is an important observation in regards to potential ecosystem consequences (e.g., decrease in SOC) during the establishment of switchgrass as a bioenergy crop in the UMRB. However, this observation was made over a single growing season with a limited sample size of soil cores. Detecting SOC is notoriously difficult because of spatial heterogeneity (Conant & Paustian, 2002), and will require further investigation with a specific sampling design before conclusions can be drawn about a potential crop-induced mechanism for SOC loss.
F I G U R E 7 Phylogenetic trees of (a) bacteria and archaea (16S)

| Microbial response to fertilization
While overall diversity was minimally affected by treatment in either the bacterial/archaeal or fungal communities, the response ratio analysis across the phylogenetic tree revealed finer-scale responses within each community. Specifically, as shown in Figure 7, the fungal community was polarized with the majority of the identified ZOTUs falling under strong positive or strong negative response classifications while bacterial/archaeal responses were centered around neutral. The lack of response in the bacterial community contradicts several long-term studies that reported opposite findings (Ai et al., 2018;Marschner et al., 2003), but this is potentially due to the lack of soil chemical changes since bacteria tend to be more sensitive to changes in abiotic factors such as pH (Strickland & Rousk, 2010). We identified several ZOTUs driving small differences between the urea and CBF:Urea treatments which primarily belonged to the Planctomycetes, Actinobacteria, and Proteobacteria phyla. While different Planctomycetes ZOTUs were enriched in either the urea or CBF:Urea treatments, the Actinobacteria and Proteobacteria ZOTUs were enriched only in the urea treatments. Long-term N fertilization has been linked to increases in the relative abundances of Actinobacteria and Proteobacteria, which has been explained by their generally copiotrophic life strategies (Dai et al., 2018;Horton et al., 2017). In addition, members of the Actinobacteria, Proteobacteria, and Planctomycetes have all been identified as potential indicators of changes in pH and C:N ratios (Hermans et al., 2017) which may suggest more microscale differences in soil chemical conditions than we were unable to detect. Furthermore, the high level of heterogeneity within soil and the legacy effects of synthetic fertilization applied to these experimental field plots may have contributed to a lack of broad community changes (Fan et al., 2019). Alternatively, since strong effects of synthetic fertilization on the microbial community tend to occur at high levels of synthetic fertilizer application (>200 kg N ha −1 ; Fierer et al., 2012;Treseder, 2008), the moderate levels used in this study (100 kg N ha −1 ) may not be sufficient to lead to significant shifts.
In the fungal community, the greater abundance of arbuscular mycorrhizal fungi (AMF) species in urea plots driving the difference between the urea and CBF:Urea treatments was surprising as the negative impact of fertilization on AMF has been reported across multiple studies and systems Ding et al., 2017;Egerton-Warburton et al., 2007;Williams et al., 2017;Zhu et al., 2018). Glomus sp. in particular have been reported to be sensitive to increasingly competitive environments such as those created by increased fertilization (Liu et al., 2015).
However, the legacy of synthetic fertilization in these experimental plots may have conditioned these members of the community so that they were less able to utilize the organic source of N from the CBF. The hyphae of AMF are known to preferentially uptake NH 4 + and also have a limited enzymatic repertoire which results in a dependence on other saprobic members of the microbial community to mineralize complex organic matter (Jansa et al., 2019). We hypothesize that a decreased ability of the fungal community to utilize the N derived from CBF biomass may have contributed to the differences in late-season NH 4 + we observed. Long-term studies comparing organic versus synthetic fertilization have reported maintenance of or increases in fungal diversity in organic systems (Kamaa et al., 2011;Kazeeroni & Al-Sadi, 2016;Liu et al., 2020), with some exceptions attributed to indirect effects of plant biomass (Piazza et al., 2019). It is thus possible that there may be a lag period for some community members during the transition from previously synthetic to organic fertilization.
In both bacterial/archaeal and fungal communities, a weak phylogenetic signal was detected, indicating that closely related species tended to respond to fertilization more similarly to each other than to randomly selected species across the phylogenetic tree. For example, although we were unable to statistically detect any specific clades that responded strongly to either fertilization treatment, it appeared that certain clades within the Firmicutes, Bacteroidetes, and Candidatus phyla responded coherently in a positive manner, while a number of archaea were negatively affected by fertilization of either type, warranting further investigation into the ecosystem functions performed by these groups. The relatively short length of this study (one growing season) and the possible effects of "relic" DNA in the soil environment on diversity measures emphasize the need for longer-term investigation of the effects of CBF on the soil microbial community (Carini et al., 2016). In addition, the contribution of the differences in abundances of individual ZOTUs to the overall community dissimilarities between treatments was exceedingly small. This suggests that even when supplementing CBF at a half-dose of synthetic fertilizer, the overriding effects of N addition on the soil microbial community are likely driven by the urea fertilization.

| CONCLUSIONS
Supplementing 50% of urea fertilizer with CBF in an outdoor field study maintained statistically similar barley biomass and increased late-season NDVI of switchgrass, indicating potential for CBF efficacy in the UMRB region. An increase in late-season N availability and evidence for mineralization of the CBF from July to September supported the potential benefits of the slower mineralization rate of CBF compared to urea. Supplementing urea fertilizer with CBF did not significantly alter bacterial or fungal composition or diversity over the growing season; however, in-depth analysis of responses across the phylogenetic tree revealed the communities had fine-scale dynamic responses. Regardless of fertilizer type, fungi had stronger responses to fertilization than bacteria with polarized strong positive and negative responses, although phylogenetic signal was detected in both kingdoms suggesting cohesive responsiveness within specific clades. In addition, we identified specific ZOTUs that had strong responses to either fertilization treatment which suggests possible avenues of further ecosystem function investigations. Importantly, the neutral response of the microbial community to biofertilizer supplementation may be construed as a lack of negative impact of the biofertilizer supplement on the soil microbial community, although longer-term studies need to be undertaken to determine if this is the case. Taken together, the results from this study are promising for the future implementation of perennial cropping with biofertilization in dryland agricultural systems.