Biological nitrogen ﬁ xation by soybean ( Glycine max [L.] Merr.), a novel, high protein crop in Scotland, requires inoculation with non-native bradyrhizobia

It is currently not recommended to grow soybean ( Glycine max [L.] Merr.) further than 54° North, but climate change and the development of new high latitude-adapted varieties raises the possibility that it could be introduced into Scotland as a novel high protein crop deriving most of its nitrogen (N) requirements through biological N ﬁ xation (BNF). This was evaluated via ﬁ eld trials in 2017 and 2018 near Dundee (56.48°N). As there are no native soybean-nodulating bacteria (SNB) in UK soils, soybean requires inoculation to exploit its BNF potential. In 2017, three commercial inoculants containing elite Bradyrhizobium strains signi ﬁ cantly increased plant biomass in plot trials with a soybean 000 maturity group variety (ES Comandor). Rhizobia were isolated from the nodules and identi ﬁ ed as the original inoculant species, B. diazoef ﬁ ciens and B. japonicum .


Introduction
Soybean (Glycine max (L.) Merr.) or soya is globally the most economically important grain legume, occupying 125 million ha and with a global production of 353 million tons per annum (data from 2020 www.fao.org,Herridge et al., 2022;Karges et al., 2022).The leading producers of soybean are the USA, Brazil, and Argentina, but with significant contributions also from China and India, although China is also a net importer (Alves et al., 2003;Herridge et al., 2022;Karges et al., 2022).The principal uses of soybean are for human consumption (oil and protein), and as a high-protein animal feed for the livestock industry.Two million tons are imported into the UK annually for animal feed alone (https://akcagric.co.uk/blog/soybean-uk-real-option; https:// www.gov.uk/government/statistics/united-kingdom-food-securityreport-2021/united-kingdom-food-security-report-2021-theme-1global-food-availability) and >3.8 m tons are consumed annually in the UK in total (Coleman et al., 2021).This represents a significant economic and environmental cost in terms of transatlantic transportation from the Americas, and the destruction of native South American vegetation to make space for its cultivation.In the UK, home-grown pulses such as fava bean have been suggested as alternatives to imported soybean, but their uptake by the animal feed industry is currently low, in part owing to their less favorable dietary attributes e.g., their lack of some essential sulphurcontaining amino acids like methionine compared to soybean, and their relatively high tannin (non-nutritional, or 'antifeedant') content (Rahate et al., 2021).There are thus a range of incentives for European countries, including the UK, to reduce soybean imports and meet their vegetable-protein demands by increasing the area of soybean grown, which in the UK currently stands at c. 10,000 ha or only 0.1% of the utilized arable land (Coleman et al., 2021).Unlike genetically modified soybean being imported from the Americas, soybean grown in Europe and the UK would be non-GM, and thus it commands a higher price and increases farmer gross profit margins (Karges et al., 2022).
Like most legumes, soybean can establish a root-nodulating symbiosis with Nitrogen (N)-fixing soil bacteria collectively termed rhizobia (Alves et al., 2003;Sprent et al., 2017).This symbiosis is capable of providing the plant with all its N-requirements through Biological N Fixation (BNF); the high BNF capacity of soybean is exploited in all the major producer countries, and has contributed to making it such a highly profitable crop by massively reducing (or even eliminating) the need for costly N-fertilizer applications (Alves et al., 2003;Schipanski et al., 2010;Leggett et al., 2017;Ciampitti and Salvagiotti, 2018;Peoples et al., 2021;Herridge et al., 2022).Soybean is generally nodulated by species of Bradyrhizobium, especially B. japonicum, B. diazoefficiens, B. yuanmingense, and B. elkanii, but in more alkaline soils by Sinorhizobium strains (Alves et al., 2003;Suzuki et al., 2008;Tian et al., 2012;Siqueira et al., 2014;Temprano-Vera et al., 2018;Zilli et al., 2021).If these symbionts are well matched with their host plant genotypes, under optimum growing conditions the soybean-Bradyrhizobium symbiosis can fix >300 kg N ha -1 year -1 (Alves et al., 2003;Hungria et al., 2005;Peoples et al., 2021), which equates to an annual global total of 25 Tg N (Herridge et al., 2022).However, there are reportedly no soybean-nodulating bacteria (SNB) in soils outside the native soybean range of Japan and China, which necessitates its inoculation in all the major areas where the crop has been introduced, including the USA, Australia, and South America.The inoculation of soybean exploits rhizobia that have most likely originated from China and/or Japan (Alves et al., 2003;Kaneko et al., 2011), but which have since become naturalized in the soybean-cropped soils of the USA and South America several decades after their introduction (Sadowsky et al., 1987;Ferreira and Hungria, 2002;Siqueira et al., 2014).Indeed, most currently available commercial soybean inoculants are based on just two strains (Siqueira et al., 2014) that are actually variants of either B. diazoefficiens SEMIA5080 (USDA122 = CB1809) or B. japonicum SEMIA 5079 (related to the USDA 123 serocluster that has become naturalized throughout the soybean-growing areas in the USA; (Cregan et al., 1989)), and are likely to be descendants of strains originally brought to the USA from Japan in the early days of soybean cropping (Van Berkum and Fuhrmann, 2000;Cytryn et al., 2008).Both strains are widely used in Brazil and Argentina, often together as a dual inoculant (Siqueira et al., 2014); in addition, B. elkanii is also sometimes used as an inoculant, solely or in combination with one of the two commercial strains (Zilli et al., 2021).
Soybean is a relatively new crop in Europe, subsequently there are no effective native soybean-nodulating bacteria (SNBs) in European soils (Madrzak et al., 1995;Albareda et al., 2009).There are relatively few studies in Europe on the yield response of soybean to inoculation with SNBs, but most field trials (although not all, e.g.Kühling et al., 2018) suggest that inoculation is generally very beneficial, significantly increasing grain and/or biomass compared to uninoculated crops (Albareda et al., 2009;Schweiger et al., 2012;Zimmer et al., 2016;Pannecoucque et al., 2018;Omari et al., 2022).In some cases, yield (and other crop benefits) in crops inoculated with SNBs are at least as good, if not better, than those receiving an N-fertilization rate as high as 200 kg N ha -1 (Albareda et al., 2009;Pannecoucque et al., 2022).In terms of BNF by the introduced soybean as measured directly using methods like N-balance, 15 N isotope dilution or 15 N natural abundance (rather than inferred from grain or biomass yields), there are even fewer studies, with most data derived from plots rather than from field-scale studies.Nevertheless, current data suggest that soybean can derive a high proportion of its N-requirements through BNF i.e., up to 60% N derived from the air (or %Ndfa; Unkovich et al., 2008) in various European countries, including Switzerland (Oberson et al., 2007), Austria (Schweiger et al., 2012), Germany (Zimmer et al., 2016), and Belgium (Pannecoucque et al., 2022).These studies demonstrate that soybean BNF in Europe, as in other parts of the World (Alves et al., 2003), is dependent on effective inoculants and low soil mineral N; whereas, with poorly-performing inoculants and/or unsuitable edaphic factors (such as high N or low P), the %Ndfa can be much lower, or even negligible.
Given the likelihood that soybean cropping will extend into Scotland in the near-to medium-future (Coleman et al., 2021), and that research into improving the cold tolerance of the symbiosis is considered to be essential for the widespread extension of soybean into Europe (Zimmer et al., 2016;Pannecoucque et al., 2018;Karges et al., 2022), this study assesses the current feasibility of growing early maturity group ("000") soybean varieties in an arable farm in East Scotland, and quantifies the resulting yield in terms of either biomass (for silage) or grain.There appears to be an absolute requirement for inoculation in all European soils that have grown soybean, and as there are mixed reports about the efficacy of commercially available products in Europe (Zimmer et al., 2016), we aim to determine the potential for soybean to be inoculated with a range of non-native Bradyrhizobium strains, and to quantify the amount of N that nodulated soybean captured through BNF under Scottish conditions using the delta 15 N natural abundance technique (Unkovich et al., 2008;Maluk et al., 2022).Finally, we assess whether introduced (non-native) inoculant bradyrhizobia could persist in Scottish soils, and what the impact of their widescale introduction and subsequent naturalization might be on the native soil microbiome.

Weather, soil properties, and experimental layout
Two field trials with soybean were conducted in consecutive years (2017,2018) at two sites more than 1 km distant from each other (South America field in 2017, House field in 2018) on the experimental farm of the James Hutton Institute at Balruddery near Dundee, Scotland, UK (56.48 lat.,.Balruddery experiences a temperate and oceanic climate.The weather conditions in East Scotland were obtained from https:// cosmos.ceh.ac.uk/sites/BALRD, and are shown in Figure S1.Insufficient precipitation in 2018, verging on drought conditions from July onwards, necessitated crop irrigation.Soil chemical properties assessed in January of both years are given in Table S1.The soil in both fields was in the Balrownie series, with a land capability rating of 3.2 (Macauley scale); the dominant soil texture was a sandy silt loam (soil depth 35 cm) with 15% stone content.

Trial 1. 2017
Seeds of the "000" very early maturity group variety ES Comandor (Euralis, France) were sown on the 23/05/2017 into plots measuring 6.25 x 1.55 m and harvested on 9/10/2017.The trial was designed to test the effects of three commercial rhizobial inoculants on soybean biomass yield, and the rate of BNF, within a total of 32 soybean plots divided into eight 1.55 m columns and four 6.25 m rows, each separated by barley (Hordeum vulgare L.) guard strips (Figures 1A, S2A).Each column included four plots (replicates) of the same treatment.There were eight plots of the same treatment in total.Four treatments were applied to seeds before sowing: control (no rhizobia application), and the following commercial inoculant products applied according to the suppliers' instructions: a) LiquiFix (http://www.legumetechnology.co.uk/) b) Rizoliq Top (http://www.rizobacter.com/brasil/produtos/inoculante/rizoliq-top/) c) Euralis (http://www.euralis-semences.fr/produits/soja/inoculum.html) The LiquiFix and Rizoliq Top applications are liquid formulations, while the Euralis product is peat-based.

Trial 2. 2018
Two soybean varieties were sown on 10/05/2018 (ES Comandor, and another '000' maturity group variety, ES Navigator) into plots/large strips.Based on the results of Trial 1, only a single commercial inoculant product, Rizoliq Top, was used in Trial 2, and it was compared with uninoculated plants.In addition to harvesting the green biomass of soybean (as in Trial 1), we also assessed the potential for yielding grain.Plots with considerably larger areas (1.55 x 100 m) were established for green biomass production (to bale for silage) together with standard-sized plots (1.55 x 6.25 m) for grain production.Trial 2 thus included a total of 24 soybean plots (12 standard plots and 12 large "strip" plots) divided into 12 columns and two 6.25 m rows, each separated by barley guards (Figures 1B, S2B).

Calculation of BNF using the 15 N natural abundance technique
In Trial 1, weather conditions did not allow for the development of pods, so the plots were harvested only for their green biomass (i.e., as potential forage).Just prior to biomass harvest, three randomly selected individual crop plants (whole shoots) were sampled per plot for d 15 N and %N analyses; these were used to estimate crop BNF using the 15 N natural abundance assay according to Maluk et al. (2022) B represents the d 15 N of soybean shoots grown entirely reliant upon BNF for growth (Unkovich et al., 2008); generation of Bvalues for soybean varieties ES Comandor and ES Navigator was performed according to Unkovich et al. (2008) and Maluk et al. (2022); all plants were inoculated with Rizoliq Top and were harvested at early-to mid-podfill stage.The non-fixing reference plants were non-nodulated soybean from the uninoculated plots.
In Trial 2, higher air and soil temperatures plus the application of irrigation ensured that the soybean crop could reach full maturity, so it was harvested for dry grain.An earlier, separate sampling aimed at capturing %Ndfa at early to mid-podfill, when grain legume BNF is maximal (Unkovich et al., 2008), was conducted on 13/09/2018.For this, soybean plants with an intact root system (to inspect for presence of nodules) were sampled for 15 N natural abundance assays and %Ndfa determinations according to Maluk et al. (2022) using the B-values obtained in Trial 1 (Table S2); with non-nodulated plants from the uninoculated plots again being used as non-fixing references.At final harvest on 23/10/2018 the mature grains were collected and their total dry weight per plot and their thousand grain weight obtained by taking into account a 12% moisture content.Sub-samples of the harvested grains were subject to d 15 N and %N determinations, and these were used to calculate %Ndfa of the grains using a soybean dry seed B-value of -0.66 delta units obtained from the literature (Bergersen et al., 1992).The amount of N fixed by the crop was then calculated from estimates of %Ndfa, shoot dry biomass (or dry grain weight) and N content (%N) using equations [2] and [3]: Legume shootðor grainÞN = % N=100 Â (legume shoot or grain dry matter) ½2 Amount shootðor grainÞN fixed = % Ndfa=100 Â (legume shoot or grain N) ½3

Rhizobia isolation and identification
To identify the rhizobial strains present in the nodules at the time of harvest in Trial 1, isolation was performed on the roots of a randomly selected soybean plant from each of the 32 plots.Fresh nodules were washed with tap water to remove residual soil and stored frozen at -80°C in 1.5 mL microtubes till needed.Isolation of rhizobial symbionts was performed from all field plots according to Maluk et al. (2022).DNA extraction from the isolates, and PCR amplification of the 16S rRNA, recA, gyrB, glnII, nodA, nodD and nifH loci was performed according to Maluk et al. (2022) with the primers used and conditions followed listed in Table S3.All Sanger sequencing was performed by the James Hutton Institute Sequencing Facility.The resulting sequences were then manually edited ('trimmed') using BioEdit Sequence Alignment Editor version 7.2 (Hall, 1999).Thereafter the sequences were compared to all publicly available nucleotide sequences in the GenBank database of the National Center for Biotechnology Information (NCBI; https://www.ncbi.nlm.nigh.gov)website using the standard nucleotide basic local alignment search tool (BLAST+ version 2.13.0) optimized for highly similar sequences (Megablast) (Altschul et al., 1997).
All gene sequences were deposited in GenBank with accession numbers OQ026739 to OQ026751 for 16S rRNA, OQ032855 to OQ032867 for recA, OQ032829 to OQ032841 for gyrB, OQ032790 to OQ032802 for glnII, OQ032842 to OQ032854 for nodA, OQ032816 to OQ032828 for nodD and OQ032803 to OQ032815 for nifH.Reference sequences were also downloaded from the NCBI database with accession numbers provided in Table S4.
All phylogenetic analyses were conducted in MEGAX (Kumar et al., 2018).First, the nodA, nodD and nifH sequences were aligned using Clustal W (Sievers et al., 2011), whereafter the best-fit evolutionary models were determined.Next, gene trees were inferred using the Maximum Likelihood (ML) approach based upon the Tamura 3-parameter model with support for the nodes being determined based upon bootstrap analyses of 1000 pseudoreplicates.All positions with less than 95% site coverage were eliminated (partial deletion option).A discrete Gamma distribution was used to model evolutionary rate differences among sites (5 categories) in the phylogenies of both the nodD (+G, parameter = 0.4368) and nodA genes (+G, parameter = 0.5903).The evolutionary model incorporated to generate the nifH phylogeny included a parameter for rate variation ([+I]), that allowed for some sites to be invariable (64.42% sites).
Sequences of the 16S rRNA, recA, gyrB and glnII loci were first used to construct single loci phylogenies with best models (T92+G +I models for 16S rRNA, recA and gyrB; TN93+G+I for glnII; data not shown).Single loci phylogenies were congruent and, therefore, could be concatenated into a multi-locus alignment, which was then analyzed (as described above) to generate the species tree.The concatenated gene tree was inferred following the ML approach and incorporating the General Time Reversible model [1].A discrete Gamma distribution was used to model evolutionary rate differences among sites (5 categories (+G, parameter = 0.5667)).The rate variation model allowed for some sites to be invariable ([+I], 69.08% sites).
Determination of survival of inoculant rhizobia in the soil using soybean trap plants and a Bradyrhizobium nodZ qPCR Soil from Trial 2 was sampled from all plots on the 04/03/2019, approximately five months after the plants were harvested.Five subsamples (300 g each) were taken from the center of each plot at equidistant intervals at a depth of 0.15 m using a soil auger or trowel (the trowel was washed with water and ethanol between plots); these were mixed and pooled into a single sample per plot.The soil samples were then stored at -80°C until they were required for qPCR analysis (see below).Sub-samples (200 g) from three of the plots per treatment (inoculated or uninoculated, plus the barley guard rows) were thawed, placed into pots, and soybean cv.Navigator seeds (three per pot) were sown into them as "trap" plants.Seeds were also sown into soil from adjacent fields without any history of soybean cropping, including those from the Centre for Sustainable Cropping (Maluk et al., 2022), but also from a field in Northumberland, England (opposite Holy Island by the causeway; 55.67 lat., -1.80 long.).Soil from a non-agricultural site in Angus, Scotland (Carrot Hill; 56.55 lat., -2.88 long.),dominated by native legumes that nodulate with Bradyrhizobium i.e. gorse (Ulex europaeus L.) and broom (Cytisus scoparius L.) (De Meyer et al., 2011;Stępkowski et al., 2018), was also examined.The trap plants were grown for 6 weeks in a controlled environment (28°C, 12 h photoperiod) during which they were given only distilled water.At harvest, their dry weights were determined, and any nodules counted and weighed.Three nodules per plant were also taken for isolation of their rhizobia according to Maluk et al. (2022).The rhizobial isolates were identified by comparing their BOX-PCR patterns with strains known to be components of commercial soybean inoculants, including those used in the present study, such as B. diazoefficiens SEMIA 5080 and B. japonicum 5079 (Siqueira et al., 2014), according to Jorrin et al. (2021).Isolated DNA was used as a template to generate BOX-PCR fingerprints, using the specific BOXA1R primer (CTACGGCAAGGCGACGCTGACG) (Versalovic et al., 1994).Amplification was carried out in a 25-ml PCR reaction containing 5-10 ng of isolated DNA and 1 U of OneTaq ® DNA Polymerase (NEB #M0480, UK).BOX-PCR products were visualized on 2% agarose gels at 100 V until clear band separation.Gel images of the BOX-PCR fingerprint of each strain in the IU population were compared to classify into operational taxonomic units (OTUs).
For the qPCR analysis, nodZ primers specific to SNBs in the B. japonicum-B.diazoefficiens group were designed and tested for their specificity.BnodZF2 (TCGTCCTCGAGCAGGTTTCGGTTAA) and BnodZR2 (CGAAGCCATAAGCGCTTGCGAGT) were found to be specific to  T and Rhizobium leguminosarum TRX19 DNA.The composite soil samples were weighed and passed through a sieve with a mesh size of 10 mm.Prior to being dried overnight at 70°C to determine their moisture content, small sub-samples of the fresh soil (2 mL in volume) were flash frozen in liquid N 2 and stored at -80°C for subsequent molecular analyses.DNA was extracted from ~0.25 g soil samples using the DNeasy PowerSoil Kit (Qiagen, Germany) following the manufacturer's instructions, using a 2 min bead-beating step (Retsch MM300, Germany) at a frequency of 30 beats s -1 .A 194 base pair (bp) fragment of the mutated 16S rRNA gene from Escherichia coli, which was used as an artificial reference spike to account for DNA losses during DNA extraction from the soil and/ or partial PCR inhibition (Daniell et al., 2012), was added to each soil sample at a concentration of 10 9 copies per sample prior to DNA extraction.
The wild-type calibrator controls were generated by PCR amplification of nodZ gene fragments from genomic DNA of B. diazoefficiens SEMIA 5080.Amplifications of the artificial spike, the nodZ calibrator controls and the nodZ gene were performed in triplicate using the LightCycler ® 480 SYBR Green I Master Kit (Roche # 04707516001, UK) following the manufacturer's instructions and using the set of primers of the amplified region of each reaction i.e.Mut-F+Mut-R for the spike (Daniell et al., 2012), and B. japonicum nodZ.The amplification conditions for the spike were: an initial denaturation at 95°C for 15 min was followed by 42 cycles of 94°C for 20 s, 58°C for 30 s, 72°C for 30 s and a single acquisition step at 81°C for 5 s.The nodZ amplification conditions were as follows: an initial denaturation at 95°C for 15 min, followed by 45 cycles of 95°C for 10 s, 69°C for 10 s, 72°C for 5 s and a single acquisition step at 85°C for 5 s.Melt curve analysis was performed between 55°C and 95°C.Samples from fields that have never been cropped with soybean served as negative controls.

Impact of soybean rhizobial inoculants on the soil microbiome
DNA was extracted from the same soybean plot soil samples as used for the trapping experiments (0.25 g each).DNA was extracted using the DNeasy PowerSoil Pro Kit (Qiagen, Germany), according to the manufacturer's instructions, the resulting DNA extracts were used to create 16S rRNA gene amplicon libraries using the primers 515F-Y GTGYCAGCMGCCGCGGTAA (Parada et al., 2016) and 806R GGACTACNVGGGTWTCTAAT (Walters et al., 2016).Amplicon libraries were created by following the standard Illumina 16S rRNA amplicon sequencing protocol.PCR cycles were as follows: initial denaturation at 95°C for 3 min, 25 cycles of denaturation at 95°C for 20 s, annealing at 55°C for 10 s, and elongation at 72°C for 20 s, with a final elongation step of 72°C for 5 min.
Samples were indexed using the Nextera XT Index kit (Illumina # FC-131-1024, USA), following the manufacturer's instructions; whereafter they were pooled in equimolar concentrations, their quality checked and sequenced (2 x 300 bp paired-end sequencing) using the Illumina MiSeq platform run at the James Hutton Institute in Dundee, UK.

Statistical analyses
Plant growth was monitored and compared to controls (uninoculated soybean plants).One-way design analysis of variance (ANOVA) and Tukey's 95% confidence intervals test in GenStat for Windows 21st Edition (VSN International Ltd., UK) were used to analyse the growth parameters, BNF and biomass quality in the 2017 field trial.
A total of 4 451 636 sequences were obtained with an average of 47, 867 sequences per sample.Bioinformatics on the sequence data was performed using QIIME 2 v2017.4(Bolyen et al., 2019).The DADA2 plugin was used to pair reads, quality control sequences and determine amplicon sequence variants (ASV).A total of 21,554 unique ASV's were assigned by DADA2 and on average after quality control and chimera removal there were 47, 576 sequences per sample.The q2-feature classifier, a naïve-Bayes classifier, and the Greengenes database were used to assign taxonomy.Samples with under 10000 sequences were assumed to have failed and were removed from the analysis.This resulted in 6 samples being removed.Relative abundances were calculated on the remaining data and principal-coordinate analysis (PCoA) with Bray-Curtis dissimilarity matrices was used to visualize b diversity of both the microbial community as a whole and of any ASVs assigned a genus of Bradyrhizobium.Differences between treatments were tested using permutational multivariate analysis of variance (PERMANOVA) within the Vegan package and R statistical environment (R version 4.1.3).

Nodulation and BNF by soybean grown in Scotland
The effect of inoculation on soybean growth was clearly visible in both Trial 1 and Trial 2 with greenish-yellow shoots of the uninoculated plants compared to the dark green foliage of the plants in the inoculated plots (Figure 1).
In Trial 1 (2017) all inoculated plants examined were nodulated regardless of the inoculant used, e.g., Rizoliq Top (Figure 1C).Most of the uninoculated plants had no nodules, but the few that did were sparsely nodulated with small, white nodules, and had the yellow foliage.Inoculation increased plant shoot biomass, albeit nonsignificantly, except in the case of Rizoliq Top (Figure 2A), and this was more pronounced (and significant for all inoculants) in terms of %N (Figure 2B) and total shoot N (Figure 2C).When adjusted to field scale, the plants inoculated with Rizoliq Top produced the highest estimated biomass (10,909 kg ha -1 ; although the differences are not statistically significant; Table S5).These plants also accumulated the most N (>350 kg N ha -1 yr -1 ), in contrast to the uninoculated plants which accumulated only slightly in excess of 100 kg N ha -1 yr -1 (Table S5).In terms of BNF, the % Ndfa values for soybean cv.ES Comandor ranged from 67% for Rizoliq Top-inoculated plants to 76% for those inoculated with LiquiFix, with the Euralis-inoculated plants intermediate at 73% although the values for the three inoculants did not statistically differ (Table S5).Application of the %Ndfa values resulted in estimates of total crop N of cv.ES Comandor obtained through BNF at c. 187-243 kg N ha -1 yr -1 depending on the inoculant product (Table S5), with Rizoliq Top-inoculated plants fixing the most N, albeit non-significantly (Figure 2C).Rizoliq Top outperformed the other two inoculant products in terms of plant growth promotion and N-accumulation despite plants inoculated with it having the lowest %Ndfa (67%) of the three products (Table S5).
In Trial 2, the dry grain yield of cv.ES Comandor inoculated with Rizoliq Top, was nearly three times higher than the uninoculated controls, i.e., 0.998 t ha -1 compared to 0.362 t ha -1 to (Figure 3A; Table S5), whilst the yield of cv.ES Navigator was also significantly raised from 0.467 t ha -1 to 0.780 t ha -1 .Inoculated plants from both soybean varieties also had significantly higher thousand grain weights (144.85 g and 149.06 g, respectively) compared to untreated control seeds (103.84 g and 102.71 g, respectively) (Figure 3B; Table S5).The N-content of the grain was also increased by inoculation, both in terms of concentration (i.e. raised from 5.30 to 7.19 mg g -1 ; Table S5), and in terms of total values per ha (71.55 kg ha -1 for cv.ES Comandor and 54.82 kg ha -1 for cv.ES Navigator).At mid-podfill stage, the %Ndfa values for inoculated ES Comandor and ES Navigator were 65% and 67%, respectively (Table S5), whereas at final harvest, they were 73% and 70%, resulting in values for total seed N obtained through BNF of 38.10 kg N ha -1 yr -1 for ES Navigator and 52.31 kg N ha -1 yr -1 for ES Comandor (Figure 3C; Table S5).In terms of plant genotype, the Rizoliq Top-inoculated ES Comandor crop significantly outperformed the corresponding ES Navigator crop in terms of grain yield, total grain N, and fixed grain N (Figure 3).

Identification of rhizobia nodulating the inoculated soybean and survival of inoculant SNB in Scottish soils
Sequence comparisons of the amplified fragments for the 16S rRNA, recA, gyrB and glnII indicated that all 13 rhizobial isolates from the sampled nodules in Trial 1 belong to the genus Bradyrhizobium.Twelve of the thirteen isolates were most closely related to the strains of B. diazoefficiens (including the type strain, USDA 110 T ) according to the concatenated 16S rRNA, recA, gyrB and glnII ML phylogeny (Figure 4), as well as individual gene phylogenies (data not shown), with this group receiving 100% bootstrap support in the concatenated phylogeny (Figure 4).This group included root nodule isolates from plots inoculated with Rizoliq Top (S7, S29, S46, S56, S61), Euralis (S21, S60), and LiquiFix (S38, S57), but also from nodulated plants in the uninoculated control plots (S34, S63, S64).A single isolate from a plot (S11) that was inoculated with LiquiFix was most closely related to strains of B. japonicum; having 100% bootstrap support for its placement with the inoculant strain B. japonicum SEMIA 5079 (Figure 4).
Sequences of the symbiosis-essential genes, nodA, nodD and nifH for most strains from this study, for the most part, formed part of the highly supported B. japonicum-B.diazoefficiens group (Figures S3-S5).This group had 100% bootstrap support across all three phylogenies and included all three strains of B. diazoefficiens as well as all three B. japonicum strains together with the type strains for B. liaoningense, B. huanghuaihaiense, B.  S3, S5); however the nodD sequence for this strain was grouped together with the sequences from the remainder of the strains recovered in this study in the nodD phylogeny (100% support; Figure S4).
Soybean plants that were sown into the soil from Trial 2 that was sampled in spring 2019 were all nodulated with approximately 30 nodules per plant regardless of whether the soil they were grown in was from plots that were inoculated in 2018 (plots 2, 7, 9) or left uninoculated (plots 14, 16, 22).Soybean plants that were sown into guard row soils (G1, G2 and G3), soils from adjacent fields that had never been sown with soybean (2011, 2013and 2016;Maluk et al., 2022), as well as soil from two other sites in the UK that had also not been cropped with soybean, were all non-nodulated.Genomic DNA of the rhizobial isolates from these nodulated plants was compared to that of standard soybean inoculant strains in a BOX-PCR analysis.The BOX-PCR profiles of all the isolates were identical to that of B. diazoefficiens SEMIA 5080 (Figure 5).A RT-qPCR using B. japonicum-B.diazoefficiens nodZ primers on extracts of the same soils used in the nodulation trials confirmed that they contained this gene which is specific to SNBs (Figure S6).Estimates of the density of SNBs in the soils derived from nodZ copy numbers suggested that the soils sown with inoculated soybean harboured more than 10 4 SNBs g -1 soil dry weight, while those from the plots sown with untreated soybean seeds harboured 10 2 -10 3 SNBs g -1 soil dry weight (Figure 6); no SNBs could be detected in the guard row soils nor in soils from fields that had never been sown with soybean.

Impact of inoculant SNB on the native soil microbiome
In total, across all samples 21,554 unique Amplicon Sequence Variants (ASVs) were detected.The most abundant phylum was Proteobacteria representing approximately 21% of all sequences followed by Firmicutes (18.5%),Actinobacteriota (13.8%),Acidobacteriota (13.5%) and Verrucomicrobiota (7%) (Figure 7A).Microbial community structure varied between soils planted with soybean and bulk soil samples from a nearby non-soybean-cropped field (p< 0.001) (Figure 7B).There was, however, no difference in community structure between soils planted with the two different soybean varieties (Figures 7A, B).In total there were 6 ASVs with similarity to Bradyrhizobium and in bulk soil this genus comprised 1.4% of all sequences.Soils planted with the variety Navigator which were not treated with Rizoliq Top and soils planted with the variety Comandor which were treated with Rizoliq Top both had Bradyrhizobium make up 1.4% of the total sequences.In contrast, soils planted with the variety Navigator which were treated with Rizoliq Top and soils planted with the variety Comandor which were not treated with Rizoliq Top showed Bradyrhizobium making up 1.6% of the total sequences in these treatments (Figure 8A).There was no significant difference in the relative abundance of Bradyrhizobium between any of the treatments or between bulk soil and soils planted with soybean (Figure 8B).

Discussion
Growth, grain yield, and BNF by soybean in Scotland We have demonstrated that over two cropping seasons, soybean can grow under Scottish conditions, and that in unusually warm years (e.g., 2018) it can even yield dry grain.In the UK, current recommendations (https://www.pgro.org/soybean/)are that soybean should not be grown further than 54°North owing to its minimum average growing temperature of 22°C (Karges et al., 2022).However, increased climate warming has the potential to Concatenated phylogeny of the 16S rRNA, recA, gyrB and glnII loci based on 45 nucleotide sequences (3007 bp in the final dataset) from rhizobia isolated from soybean nodules sampled in Year 1 (2017) in comparison with already-known soybean strains and with Bradyrhizobium type strains.The tree is unrooted.The type strains are indicated by a 'T' at the end of each strain code.The trees are drawn to scale, with branch lengths measured in the number of substitutions per site.
extend soybean cropping into Northern areas of the UK, and our data suggest that this northward trend may already be underway, with drier conditions making it possible to grow soybean with harvestable dry grain.Current soybean varieties used in Europe have been bred for early maturity so that they can be harvested before the onset of cooler temperatures in early autumn (Karges et al., 2022), but it is still likely that grains will not mature in most growing seasons in the north of the UK.Therefore, growers in such marginal growing areas could consider selling green grains for edamame beans.Alternatively, if as occurred with our 2017 trial, there was no harvestable grain, but there was considerable production of foliage (estimated from plots as > 10 t ha -1 in the case of plants inoculated with Rizoliq Top), then harvesting for whole-crop silage could still be a profitable option (Vargas-Bello-Peŕez et al., 2008).
We have demonstrated that soybean can be grown in Scotland without any added N-fertiliser by largely exploiting its ability to utilise BNF from its nodulating symbiosis with Bradyrhizobium.Nodulated soybean plants inoculated with Bradyrhizobium were visibly healthier (greener) than their uninoculated counterparts in both trial years, indicating that only the inoculated plants had sufficient N to sustain the production of chlorophyll.The few uninoculated plants that were sparsely nodulated with small, white nodules had yellow foliage, indicating severe N-deprivation  Quantification of Bradyrhizobium inoculant strains that survived in soils over winter 2018/2019 using a qRT-PCR method with soybean nodulating bacteria (SNB)-specific nodZ primers.Each pair of boxplots with the same letter did not show significant differences (p-value<0.05) between each other according to the results of Tukey's HSD.The nodZ gene was amplified from soil samples inoculated with Rizoliq TOP.There was no amplification of nodZ in guard row soils (G1, G2 and G3), nor in soils from adjacent fields (2011, 2013 and 2016) and from sites in other parts of UK that had never been sown with soybean: Northumberland (a field by the causeway to Holy Island, NE England, UK) and Angus (Carrot Hill 5 km NE of Dundee, Scotland, UK). and suggesting that these nodules had contributed little or no fixed N by the time they were harvested.In the case of cv.ES Comandor, inoculated plants produced up to 36% more dry matter in 2017, and 200% more dry grain in 2018 than the uninoculated plants.The % Ndfa values derived from the 15 N NA technique in both trial years were in the range of 65 -76% regardless of plant variety and inoculant used, and are slightly higher than average global values for soybean (Salvagiotti et al., 2008;Ciampitti and Salvagiotti, 2018), but similar to %Ndfa values obtained from other European studies (Oberson et al., 2007;Schweiger et al., 2012;Zimmer et al., 2016;Pannecoucque et al., 2022).The remaining 24 -35% N in the plants was most likely due to direct root uptake.In Trial 1, plant available soil N at depths down to 90 cm was 162.4 kg ha -1 (Table S1), which would have been sufficient to "top up" the N-requirements of the inoculated/nodulated plants but was presumably insufficient to provide the uninoculated plants with enough N for healthy growth (possibly because they could not access much of the deeper N-pool below 60 cm).In Trial 2, soil available N was much lower (34.72 kg ha -1 ), but would still have been sufficient to help account for the total amount of N in the grains at final harvest of the inoculated plants, and would also explain the source of all the (greatly reduced) N in the grain of the uninoculated plants that had little or no BNF contributions.However, we did not assess total plant BNF at any earlier stages than final harvest in Trial 2, so it is likely that the fixed N in the grain represents a fraction of the total N fixed by the plant during its lifetime, with much of the N being lost to the soil and atmosphere via root exudates and the shedding of leaves (Ciampitti and Salvagiotti, 2018).
Persistence of SNB in Scottish soils and their impact on the soil microbiome The SNB isolated from nodules sampled in 2017 were almost all identified as belonging to B. diazoefficiens, with only a single B. japonicum isolate from a plot inoculated with LiquiFix.We also observed that only strains with BOX-PCR profiles almost identical to B. diazoefficiens SEMIA5080 could be isolated from the nodules on soybean "trap" plants sown into the soils sampled in spring 2019 from the plots cropped with soybean in 2018.Collectively, these data suggest that B. diazoefficiens can persist in Scottish soils over winter, i.e. for at least six months after the soybean had been harvested.This is in acordance with studies from France (Obaton et al., 2002), Poland (Narozṅa et al., 2015) and Germany (Yuan et al., 2020;Halwani et al., 2021) demonstrating that bradyrhizobial SNBs persist in Northern European soils, but is the first observation of this phenomenon in the wetter and colder climate of a country like Scotland.The persistence of SNB in both the inoculated and uninoculated soils in the present study could be linked to the presence of host plant residue in the soil (Moroenyane et al., 2020).On the other hand, the apparent persistence of specifically B. diazoefficiens in Scottish soils is surprising given that most currently available commercial soybean inoculants contain a mixture of B. diazoefficiens and B. japonicum strains (Siqueira et al., 2014), and of these, the B. japonicum strains (especially those in the highly competitive USDA123 sero group) are considered to be the more persistent SNB in important soybeancropping regions, such as the USA (Sadowsky et al., 1987;Cregan et al., 1989) and Brazil (Ferreira and Hungria, 2002;Mendes et al., 2004;Batista et al., 2007).The soil type, climate, and the genotype of the inoculant strains is likely to be of importance as has been observed in other locations like Spain (Albareda et al., 2009), Brazil (Ferreira andHungria, 2002, Silva Batista et al., 2007), Poland (Narozṅa et al., 2015), and the USA (Sadowsky et al., 1987), so it is possible that the soils in our trials have facilitated the survival of B. diazoefficiens over that of B. japonicum.Clearly, further research is needed to determine which plant, soil, and climatic factors assist in the persistence of particular Bradyrhizobium genotypes (Albareda et al., 2009;Moroenyane et al., 2020).
The trap plants sown in soil sampled in March 2019 provided strong evidence that SNB were still present in all the Trial 2 plots, regardless of treatment.However, the trap plants could not provide information about the number of SNB persisting in the Trial 2 plots over winter.Therefore, a qRT-PCR method using nodZ primers specific to SNB was employed to obtain data.This methodology, based on the counting of gene copy numbers, has been used to quantify soil populations of other legume symbionts, such as Rhizobium leguminosarum (MacDonald et al., 2011;Maluk et al., 2022) and has advantages over the "traditional" Most Probable Number (MPN) method in that it is faster, more accurate, and cheaper (Furseth et al., 2010).Although a similar approach was used to enumerate SNB in Midwestern USA soils (Furseth et al., 2010), the present study is the first to attempt to quantify SNBs in European soils using a qRT-PCR method.This method does not distinguish between SNB types or species as it counts copy numbers of the nodZ gene, a host-specific nodulation gene essential for successful recognition between soybean and Bradyrhizobium strains resulting in nodule formation.This gene is very similar, if not identical, in many Bradyrhizobium SNBs, such as B. japonicum and B. diazoefficiens (Furseth et al., 2010), but the isolation of just B. diazoefficiens from the trap plant nodules (per the BOX-PCR results) suggests that this was the predominant genotype in the soil.In the soil from the inoculated plots, nodZ was present in quite high copy numbers i.e. 10 4 g -1 dry soil.These SNB populations were similar to other studies of post-SNBinoculated European soils (albeit using the MPN or fluorescence methods); although importantly, those used soils that had much longer periods (several years) of soybean cropping (Crozat et al., 1982;Madrzak et al., 1995;Albareda et al., 2009).
For the soils sampled after the 2018 trial, SNB were detected in the uninoculated plots, albeit at numbers 10-to 100-fold less than in the inoculated plots.These uninoculated soils were capable of producing nodulated soybean in trap plant trials, indicating that their lower SNB populations were still sufficient to effectively nodulate soybean.This confirms data from Furseth et al. (2010) that as few as a single rhizobium cell (per g soil) needs to be present for nodulation to occur.The two most likely (and non-mutually exclusive) explanations for SNB being present in the uninoculated plots is that they were either derived from the few nodules that formed on the plants in the uninoculated plots by the time that the soybean were harvested, or they were transferred from the inoculated plots over the period up until soil sampling in March 2019 via wind, rain, surface water run-off, human activities, or farm vehicles (Narozṅa et al., 2015).The presence of native SNB is extremely unlikely as no nodules were formed in any of the non-soybean-cropped soils that were tested, including the barley guard rows between the plots.On the other hand, despite several studies suggesting their absence (Madrzak et al., 1995;Narozṅa et al., 2015), the potential occurrence of native SNB in Europe has been hinted at by a Belgian "citizen science" study in Flanders (Van Dingenen et al., 2022) in which a number of isolates were obtained from nodules that developed on soybean plants sown into garden soils, although their efficacy and practical application remains to be demonstrated.There is little doubt, however, that inoculant SNB can over several years become naturalized in non-native soils (Sadowsky et al., 1987), and in some environments the possibility of horizontal gene transfer (HGT) of symbiosis-related genes from inoculant SNB to native strains has been demonstrated, such as in central Brazil (Silva Batista et al., 2007).
The soil microbiome analysis suggested that the impact of the introduction of the soybean plant itself was more important than that of the SNB inoculation, which is not unusual when a novel crop is introduced into a soil (Zgadzaj et al., 2016).Most studies on the soybean microbiome (both of the plant itself and its rhizosphere) and the impact of SNB inoculation, come from China where the soils harbour high populations of native SNB (e.g., Gao et al., 2022) and so are not directly comparable to the SNB-naïve Scottish soils.For example, Zhong et al. (2019) showed that the abundance and diversity of Bradyrhizobium increased after SNB inoculation, but also resulted in an increase in the abundance of other potentially beneficial bacteria.In the present study, the occurrence of soybean increased the diversity of bradyrhizobia compared to bulk soil with no soybean; this supports the nodZ qRT-PCR analysis which demonstrated that there were SNB present in both the inoculated and uninoculated plots.Therefore, it is likely that at least some of the increased bradyrhizobial diversity in the soybean-cropped plots was due to the effect of adding SNB.The other sources of bradyrhizobia in Scottish soils are likely to be background populations of non-symbiotic bacteria (VanInsberghe et al., 2015) together with symbiotic bradyrhizobia that nodulate native Genistoid legumes like gorse and broom that do not harbour soybean-compatible nod genes (De Meyer et al., 2011;Stępkowski et al., 2018).

Conclusions
Soybean can be grown in Scotland for either biomass (for forage, silage etc.) or, weather-permitting, for dry grain.The BNF ability of nodulated soybean can be exploited to provide most of the plant's N-requirements, but it requires inoculation with non-native Bradyrhizobium strains, as failure to inoculate results in significant reductions in both yield and grain quality.The expansion of soybean cropping into Scotland could be enhanced by the introduction of much earlier-maturing varieties that are also more cold-tolerant (Coleman et al., 2021;Karges et al., 2022).In common with other northern European countries, inoculant bradyrhizobial strains can persist in Scottish soils (at least over a single winter as we have shown here), which suggests that if inoculated soybean is increasingly and more frequently cropped then SNB populations could build up in the soil to the point where annual inoculation becomes unnecessary (Obaton et al., 2002;Alves et al., 2003).However, there is an urgent need for further research and monitoring to determine the impact of the introduction of high populations of SNBs into naïve soil, such as the naturalization of highly competitive but relatively ineffective strains, the transfer via HGT of SNB symbiosis genes to native bradyrhizobia, and ecosystem perturbations, such as reductions in diversity and resilience of the native soil microbiota caused by (for example) increased competition and antagonism (Mawarda et al., 2020).

Data availability statement
The datasets presented in this study can be found in online repositories.The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.
, i.e. the proportion of soybean N derived from atmospheric N 2 (%Ndfa) was calculated by comparing the 15 N natural abundance (expressed as d 15 N or parts per thousand [‰] relative to the 15 N composition of atmospheric N 2 ) of the nodulated soybean shoot N (d 15 N legume) to the d 15 N of the non-N 2 -fixing reference plants (which are assumed to reflect the d 15 N of the plantavailable soil N [d 15 N soil]) using equation [1]: % Ndfa = 100 Â (d 15 N soil − d 15 N legume)=(d 15 N soil -B) ½1

FIGURE 1
FIGURE 1 Photos of the two field trials on the experimental farm of the James Hutton Institute at Balruddery near Dundee, Scotland, UK (56.48 lat., -3.13 long.)(A) Year 1; 2017 (B) Year 2, 2018.(C) Nodulated soybean plant inoculated with Rizoliq TOP in Year 1. (D) Nodulated soybean plant inoculated with Rizoliq TOP in Year 2.
FIGURE 2Shoot dry biomass (A; F pr. 0.054), shoot N concentration (%N) (B; F pr.<.001) and total crop N (F pr.<.001) and total crop fixed N (F pr.<.001) at mid pod-fill (C) of soybean during the 2017 growing season at the experimental farm of the James Hutton Institute at Balruddery near Dundee, Scotland, UK.Data are means ± standard error.Significant differences (p<0.05) are indicated with different letters.
FIGURE 3Grain yield (A; F pr.<.001), thousand seeds weight (%N) (B; F pr.<.001) and total crop N (F pr.<.001) and total crop fixed N (F pr.<.001) (C) of soybean seeds from the Year 2 (2018) field trial at the experimental farm of the James Hutton Institute at Balruddery near Dundee, Scotland, UK.Data are means ± standard error.Significant differences (p<0.05) are indicated with different letters.

FIGURE 6
FIGURE 6 FIGURE 7 b diversity (A), relative abundance of the top 10 most abundant phyla (B) from 16S rRNA libraries of soil sampled in spring 2019 from the site of the Year 2 (2018) field trial at the experimental farm of the James Hutton Institute at Balruddery near Dundee, Scotland, UK.

b
FIGURE 8 b diversity of Bradyrhizobium (A), relative abundance of Amplicon Sequence Variants (ASVs) identified as Bradyrhizobium (B) from 16S rRNA libraries of soil sampled in spring 2019 from the site of the Year 2 (2018) field trial at the experimental farm of the James Hutton Institute at Balruddery near Dundee, Scotland, UK.