Microbiome structure and response to watering in rhizosphere of Nitrosalsola vermiculata and surrounding bulk soil

The plant rhizosphere microbiomes were thought to help the plant stands adverse condition. The study aims at deciphering signatures of rhizosphere soil microbiomes of the medicinal plant Nitrosalsola vermiculata and those of the surrounding bulk soil as well as to detect influence of watering in restructuring soil microbes that can improve the plant’s ability to tolerate drought stress. Amplicon sequencing of partial 16S rRNA gene indicated that alpha diversity indices are higher in rhizosphere than in bulk soils, while no distinctive differences were observed due to the watering. Relative abundance of phylum Cyanobacteria and its descendent unidentified genus is the highest among phyla and genera of bulk soil. Relative abundance of phyla Euryarchaeota, Chloroflexi, Actinobacteria, Proteobacteria, Bacteroidetes, Firmicutes, Acidobacteria and Gemmatimonadetes as well as genera Bacillus, Ammoniphilus, Sphingomonas, Microvirga, Pontibacter, Adhaeribacter and Arthrobacter was significantly higher in rhizosphere soil. The latter taxa were reported to act as plant growth-promoting bacteria (PGPB) through symbiotic associations. We speculate that relative abundance and mutual dominance of these taxa in rhizosphere of N. vermiculata were due to the intensity and type of plant root exudates. Other factors include soil pH where microbes favoring high soil pH can show better growth in rhizosphere soil. Also, co-existence of phyla that promote sustainability of cohabiting phyla in the rhizosphere and have high synergism prevalence in biofilm formation can be one extra factor. Quorum sensing (QS) also mediates bacterial population density in a given environment and elicit specific plant responses. The low abundance of Cyanobacteria in rhizosphere soil can be due to the inhibitory effect of highly abundant members of Firmicutes, especially those of genus Bacillus. The latter conclusion was confirmed by the occurrence of high expression rate of comQ gene triggering QS in genus Bacillus. Highly abundant microbes whose abundance was not changed due to watering are phyla Firmicutes, Proteobacteria, Chloroflexi and Cyanobacteria and their descendent genera Bacillus, Ammoniphilus, Sphingomonas, Microvirga and unidentified genus of Cyanobacteria. We speculate that non-responsive taxa to watering were drought tolerant AcademicPres Notulae Botanicae Horti Cluj-Napoca Agrobotanici Abuauf HW et al. (2022). Not Bot Horti Agrobo 50(1):12567 2 and can help plants stand adverse conditions of water scarce. In conclusion, insights on the factors involved in shaping microbiome signatures and those eliciting differential plant responses to drought stress are raised and warrant further investigations.

geochemistry, e.g., magnitude and directionality of soil community shifts (Naylor et al., 2017). Shift magnitude is shown to be proportional to the strength and duration of drought, which will be diminished as soon as water becomes available (Xu et al., 2018). This indicates that patterns of bacterial growth are dynamic and fastresponding to adverse environmental conditions. Several reports indicate that root microbiome shifts under drought stress are in favor of phylum Actinobacteria and other Gram-positive taxa (Naylor et al., 2017;Timm et al., 2018;Xu et al., 2018). It is previously argued that Gram-positive (diderm) bacteria are more drought tolerant than Gram-negative (monoderm) bacteria because the first have thicker cell wall and harbor several strategies for drought avoidance (ex., spore formation) (Potts, 1994). We suspect that many other intermicrobial and plant-microbial relationships as well as type of soil might influence drought tolerance of plants in arid regions.
The aims of the present study were to detect microbiome structures in rhizosphere soil of N. vermiculata and surrounding bulk soil as well as to detect soil microbial shifts and restructuring due to water.

Materials and Methods
Watering experiment and soil collection Watering experiment was conducted in the northern region of Mecca, Saudi Arabia, where N. vermiculata perennial plants grow naturally in the wild (Al-Eisawi and Al-Ruzayza, 2015). A spot that received no rainfall for at least three months prior experiment was selected for soil sample collection (~10-30 cm depth) as previously described (Dai et al., 2019;Geng et al., 2018). Assigned plots (1 m 2 each) for the experiment contain single-grown similar-sized plants in addition to plots of neighboring bulk soil. Three plots were watered once in the morning (25 L dH2O/plot), then, rhizosphere and bulk soil samples were collected after 0, 24 and 48 h of watering in three replicates. Soil samples were daily collected immediately put in liquid N, transported on dry ice, and stored at -20°C as described (Hurt et al., 2001). Then, soil pH was determined in a 1:1 (wt/wt) soil-H2O slurry as described (Mclean, 1983) and samples were split for further extraction of metagenomes and metatranscriptomes.
Metagenomic DNA extraction and amplicon sequencing DNAs of the different samples were extracted using CTAB/SDS method and concentration and purity were checked by electrophoresis on 1% agarose gels. DNA concentration was adjusted to 1 ng/μL prior shipment to Novogene Co., Ltd., Singapore for deep sequencing of the partial 16S rRNA gene at region V3-V4 was done. Then, libraries were generated using NEBNext® Ultra TM DNA Library Prep kit for Illumina platform and deep sequencing resulted in the recovery of ~250 bp paired-end raw reads that were quantified via Qubit. Then, bioinformatics analysis was done as described below.

Diversity and taxonomy analyses
Raw read datasets were firstly generated for microbiome samples of surrounding bulk (So) and rhizosphere (Rh) soils of N. vermiculata were collected in three replicates at the three watering time points as following: 0 (So11-So13 or group A for bulk soil & Rh11-Rh13 or group D for rhizosphere soil), 24 (So21-So23 or group B for bulk soil & Rh21-Rh23 or group E for rhizosphere soil) and 48 h (So31-So33 or group C for bulk soil & Rh31-Rh33 or group F for rhizosphere soil). Further nomenclature includes grouping style ABCDEF representing interaction between soil type and time after watering, e.g., 0, 24 and 48 h. Raw datasets were also grouped based on soil type regardless of watering time points, e.g., bulk (group J) and rhizosphere (group K) soils comprising grouping style JK, and on time after watering, e.g., 0 (group G), 24 (group H) and 48 h (group I) regardless of soil type comprising grouping style GHI.
Then, raw datasets were merged using FLASH (V1.2.7) (Magoč and Salzberg, 2011) and high-quality clean tags (Bokulich et al., 2013) were subjected to quality filtering via QIIME software V1.7.0 (Caporaso et al., 2010). The recovered tags were compared with the reference database (SILVA database, http://www.arbsilva.de/) and chimeric sequences were removed via UCHIME algorithm (Edgar et al., 2011). Effective tags were analyzed using Uparse software V7.0.1090 (Edgar, 2013), where sequences with ≥97% similarity were considered as a single OTU. Then, QIIME software V1.7.0 in Mothur method (Altschul et al., 1990) was used in combination with SSUrRNA dataset of SILVA database for taxa annotation (threshold of 0.8-1) and for assigning taxa to different taxonomic levels, e.g., kingdom, phylum and genus (Quast et al., 2012). Ternary plots were generated to detect the 10 most dominant taxa at phylum and genus levels among group style ABC for bulk soil microbiomes, DEF for rhizosphere microbiomes and GHI for microbiomes at the three time points as described (Bulgarelli et al., 2015).

Alpha and beta diversity analyses
Alpha diversity indices of Shannon and Simpson were generated with QIIME software V1.7.0 and displayed via R software V2.15.3 (Li et al., 2013). Statistical analysis for boxplots of the two indices was formed between individual groups of grouping styles ABCDEF, GHI and JK using t-test. Rarefaction curves were 4 performed, sequencing data volume was detected for rationality, and species accumulation boxplot or curve (Specaccum) was generated as described (Lundberg et al., 2013).
Anosim test of beta diversity was primarily performed by R software (Vegan package: anosim function), then, other beta diversity indices were estimated to detect the differences between or among microbial communities using weighted and unweighted unifrac matrices of QIIME software V1.7.0. The data matrices were used in generating heatmaps of weighted and unweighted unifrac beta diversity, in measuring unweighted pair-group method with arithmetic means (UPGMA) and in measuring principal coordinate analysis (PCoA) to obtain ecological distances among samples. UPGMA method was used to construct phylogenetic trees with relative abundance by phylum (Lozupone and Knight, 2005;Lozupone et al., 2011;Lozupone et al., 2007) and PCoA was displayed by WGCNA package, stat packages and ggplot2 package in R software V2.15.3. Evolutionary phylogenetic trees for individual samples and their grouping styles at genus level (the top 100 genera in abundance) were constructed using MUSCLE V3.8.31 (Edgar, 2004). Then, phyla and genera at levels of grouping styles ABCDEF, GHI and JK were selected for further analysis based on statistical significance.
Extraction of soil total RNAs and real time PCR Soils of the three replicates of individual groups of grouping style ABCDEF were gathered. Then, total RNAs were isolated using RNA PowerSoil® Total RNA isolation kit (Mo Bio, cat. no. 12866-25) following manufacturer manual. DNA contamination was removed using RQ1 RNase-free DNase (Promega, Madison, WI, USA). Absence of detectable carryover DNA in the processed samples was confirmed via agarose gel electrophoresis after performing regular PCR of comQ partial length gene of Bacillus subtilis without prior reverse transcription amplification. Primers to amplify comQ (NC_000964.3) gene fragment (FW: 5` TTGCATGGCCTCGTTTACG 3`, & RV: 5` CTCCTTTGCTGAATCCACAATC 3`, 289 bp) of B. subtilis as well as its 16S rRNA (AB042061) house-keeping gene fragment (FW: 5` TCGCGGTTTCGCTGCCCTTT 3`, & RV: 5` AAGTCCCGCAACGAGCGCAA 3`, 177 bp) were designed using Netprimer software (http://www.premierbiosoft.com/netprimer/index.html) following standard criteria. Expression level of the metatranscript comQ was detected by real time PCR using Agilent Mx3000P qPCR System (Agilent technology, USA). Maxima TM SYBR Green/ROX real time PCR was done as described (Bahieldin et al., 2015). Similar amounts of RNAs were used for the six group samples and calculations were made to detect the expression level of comQ gene in individual samples relative to that of the house-keeping gene.

Results
Deep sequencing of 16S rRNA partial-length gene was performed for microbiota of rhizosphere soil of the perennial plant N. vermiculata and surrounding bulk soil in order to detect differential growth dynamics in the two soil types and those due to watering. The amount of water used (25 L dH2O/m 2 plot) kept the soil moist during sample collections across the three time points. We thought microflora will respond instantly to watering by changing its diversity and interactions with the intact plant.
Statistics of raw sequencing data We used Illumina MiSeq for analyzing microbiomes of rhizosphere soil of N. vemiculata versus surrounding bulk soil after watering at three time points 0, 24 and 48 h. Raw sequencing data are shown in Table 1. Microbiome grouping styles refer to the interaction between soil type and watering time point (e.g., grouping style ABCDEF), to time points after watering regardless of soil type (e.g., grouping style GHI), and to soil types regardless of watering time points (e.g., grouping style JK). Sequence length per read ranges between 409-416 bp with average raw, clean and effective tag numbers of 188225, 184389 and 146150, respectively. Mean percentages of Q20 and Q30 with sequencing error rates of ˂ 1 and 0.1% are 97.56 and 92.69%, respectively, and mean percentage of effective tags is 70.67% (Table 1). The data in Figure 1 indicate that the mean number of OTUs per sample is 2273, while numbers of total, taxon and unique tags are 146150, 138889 and 7258, respectively. Expectedly, average number of archaeal sequences across samples is as low as 672, while that of bacterial sequences is 138217 ( Figure S1). Table 1. Statistics of sequencing data of microbiomes collected from surrounding bulk (So) and rhizosphere (Rh) soils of Nitrosalsola vermiculata in three replicates after 0 (So11-So13 or group A & Rh11-Rh13 or group D, respectively), 24 (So21-So23 or group B & Rh21-Rh23 or group E, respectively) and 48 h (So31-So33 or group C & Rh31-Rh33 or group F, respectively) of watering. Samples data were also grouped based on source of microbiome regardless of time after watering, e.g., bulk (group J) and rhizosphere (group K) soils and on time after watering regardless of source of microbiome, e.g., 0 (group G), 24 (group H) and 48 h (group I) Raw PE no. represents number of original paired-end (PE) reads after sequencing; Raw tags no. represents number of tags merged from PE reads; Clean tags no. represents number of tags after filtering; Effective tags no. represents number of tags after filtering chimera and can be finally used for subsequent analysis; Base nt no. is the number of bases of the Effective Tags; Avg. len. (nt) represents average length of Effective Tags; Q20 and Q30 are the percentages of bases whose quality value in Effective tags is greater than 20 (sequencing error rate is less than 1%) and 30 (sequencing error rate is less than 0.1%); GC (%) represents GC content in Effective Tags; Effectivity (%) represents the percentage of Effective Tags in Raw PE.
In terms of bacterial description across taxonomic ranks, Figure S2 shows that annotation ended at genus level is, expectedly, richer than those of either higher or lower ranks. Table S1 refers to tagged OTUs and their taxonomic description (kingdom, class, order, family, genus and species). The total number of OTUs across samples is 4180 of which number of archaeal OTUs is 36. Number of OTUs with average number of sequencings reads over 500 is 20 of which all except one taxon belong to phyla Cyanobacteria, Proteobacteria, Fimicutes, Bacteroidetes and Actinobacteria, while the highest abundant taxa of phyla Acidobacteria (OTU31), Chloroflexi (OTU74) and Gemmatimonadetes (OTU183) have as little as ~ 250, 236 and 95 sequencing reads, respectively (Table S1). The largest average number of reads is assigned to OUT1 (~32147) of unidentified genus of Cyanobacteria, followed by OTU2 (~7416) of an unidentified bacterial phylum, OTU5 (~5726) of genus Microvirga, and OTU7 (~4712) of genus Bacillus (Table S1). Interestingly, average number of reads for the unidentified bacterial phylum in bulk soil (14574) is far higher than that in the rhizosphere (258). This indicates that bulk soil has a wealth of microbes remains to be deciphered at the phylum level. The top 10 most abundant phyla and genera in grouping styles ABCDEF and JK are shown in Figure 2, while grouping style GHI is shown in Figure 3. Aligning with our previous results, relative abundance of phylum Cyanobacteria was the highest among different phyla in grouping styles ABCDEF and JK, followed by the eight taxa of bacterial phyla Proteobacteria, Firmicutes, Actinobacteria, Bacteroidetes, Acidobacteria, Chloroflexi, Gemmatimonadetes and the unidentified bacteria, as well as the archaeal phylum Euryarchaeota. Figure 3 indicates that unidentified genus of Cyanobacteria was the highest among different genera in grouping styles ABCDEF and JK, followed by the nine genera Bacillus, Microvirga, Adhearibacter, Arthrobacter, Pontibacter, Ammoniphilus, Sphingomonas, Bacteroidetes in addition to the highly abundant unidentified bacteria.

Alpha diversity analysis
Results of alpha diversity at the individual sample level indicate that species richness and evenness are higher in rhizosphere than in bulk soils ( Figure S3). Similar trend of results is reached at the grouping style level, where Shannon and Simpson indices are higher in grouping style DEF (or group K) than in grouping style ABC (or group J) ( Figure 4). In terms of microbiomes collected at different time points of watering across soil types (e.g., grouping style GHI), no distinctive differences in both Shannon and Simpson indices are observed. Rarefaction curves show that the maximum depth permitted to retain all samples in the dataset is ~120,000 sequence reads ( Figure S4). While, species accumulation curve or Specaccum indicates continuous increase of the curve height as number of sequenced samples increases (up to 18) ( Figure S5). Towards the end of the graph, the curve is flattened indicating that microbiomes of the two soil types are practically saturated. 7 Figure 2. Relative abundance of the top 10 taxa at the phylum and genus levels within grouping styles ABCDEF and GHI of microbiomes collected from surrounding bulk (grouping style ABC) and rhizosphere (grouping style DEF) soils of Nitrosalsola vermiculata after 0 (grouping style AD or group G), 24 (grouping style BD or group H) and 48 h (grouping style CF or group I) of watering . Not Bot Horti Agrobo 50 (1) Figure S6). On the other hand, anosim boxplot between groups J and K shows an R value of 1 (P = ~0.00) indicating that inter-group differences are significantly greater than intra-group differences. Anosim is a nonparametric test displayed as boxplots to justify grouping styles and to evaluate whether variation among groups is significantly higher or lower than variation within groups. Therefore, the results of grouping style JK are justified by anosim and indicate that grouping style JK referring to type of soil has high rationality for grouping ( Figure S6), while grouping style GHI referring to the watering time points has low rationality. Figure 5 refers to heat maps of weighted and unweighted unifrac beta diversity measures at the different grouping style levels. Weighted uniFrac describes relative abundance of microbes, while unweighted uniFrac refers to the unique species in the environmental samples (Lozupone and Knight, 2005;Lozupone et al., 2007;Lozupone et al., 2011). The results of heat maps of grouping styles ABCDEF and JK indicate high beta diversity distance between microbiomes of bulk and those of rhizosphere soils in terms of weighted unifrac, while moderate beta diversity distance in terms of unweighted unifrac ( Figure 5). The heat map in grouping style GHI almost indicate no distance in weighted unifrac among groups G, H and I, while moderate distance in unweighted unifrac matrices. Then, we expect to get more enriched micrbiomes in grouping styles ABCDEF and JK than in grouping style GHI. Figure S7 refers to dendrogram trees or hierarchical clustering by phylum of unweighted pair-group method with arithmetic mean (UPGMA). 9 Figure 4. Boxplots of Shannon and Simpson alpha diversity measures referring to taxa richness and evenness, respectively, of microbiomes collected from surrounding bulk (grouping style ABC) and rhizosphere (grouping style DEF) soils of Nitrosalsola vermiculata after 0 (grouping style AD), 24 (grouping style BE) and 48 h (grouping style CF) of watering Samples data were also grouped based on time after watering regardless of soil type, e.g., 0 (G), 24 (group H) and 48 h (group I) and on soil type regardless of time after watering, e.g., bulk (group J) and rhizosphere (group K) Figure 5. Heatmaps of weighted (top record) and unweighted unifrac beta diversity measures of microbiomes collected from surrounding bulk (grouping style ABC) and rhizosphere (grouping style DEF) soils of Nitrosalsola vermiculata after 0 (grouping style AD), 24 (grouping style BE) and 48 h (grouping style CF) of watering Samples data were also grouped based on on time after watering regardless of soil type, e.g., 0 (G), 24 (group H) and 48 h (group I) and soil type regardless of time after watering, e.g., bulk (group J) and rhizosphere (group K) soils. Figure 5 indicates descrete separation among samples of bulk and rhizosphere soils in terms of weighted unifrac measure, while separation in terms of unweighted unifrac measure is low. This demonstrates that microbial relative abundance at different grouping style levels, especially JK, is feasible. UniFrac refers to the phylogenetic information of environmental samples, but when coupled with standard multivariate statistical techniques, it refers to principal coordinates analysis (PCoA). The latter maps the distances among microbial communities. The results in Figure 6 indicate complete separation among microbiomes of bulk and rhizosphere soils in grouping style JK, as the first was located in the negative position of PCoA 1 (or PC1), while the second was located in the positive position of PC1. The results for grouping style ABCDEF show complete separation of the microbiomes of rhizosphere soil in terms of the three watering time points (groups D, E and F). Bulk soil microbiome of group A referring to 0 h time point of watering shows complete separation in the negative positions of PC1 and PC2, while those of groups B and C referring to 24 and 48 time points of watering merge in the negative and positive directions of PC1 and PC2, respectively ( Figure 6).
Ternary plots describing the differences in dominant phyla in subgrouping styles ABC of bulk soil and DEF of rhizosphere soil indicate that phylum Cyanobacteria, followed by Proteobacteria, is the most dominant in bulk soil microbiomes, while Proteobacteria, followed by Actinobacteria, in rhizosphere soil ( Figure S8). At the genus level, unidentified genus of Cyanobacteria is the most dominant in bulk soil microbiomes, while genus Bacillus followed by Microverga, in rhizosphere soil ( Figure S9). In terms of ternary plots in grouping style GHI across soil types, phylum Cyanobacteria again, followed by Proteobacteria, while unidentified genus of Cyanobacteria, followed by genus Bacillus, are the most dominant taxa ( Figure S10). watering Samples data were also grouped based on time after watering regardless of soil type, e.g., 0 (group G), 24 (group H) and 48 h (group I) and on soil type regardless of time after watering, e.g., bulk (group J) and rhizosphere (group K) soils. A colored dot represents a given sample in one group, and similar colored dots refer to samples of the same group. X-axis is the first principal coordinate and Y-axis is the second. Number in brackets represents contributions of PCoAs to differences among samples.
The results of ternary plots in grouping styles ABCDEF and GHI at the phylum and genus levels support our previous findings in Figures 2 and 3. Evolutionary phylogenetic trees of representative sequences by genus (the top 100 genera in abundance) are constructed for samples ( Figure S11) and their grouping styles ABCDEF ( Figure S12), GHI ( Figure S13) and JK ( Figure S14). The results of phylogenetic trees for samples and different grouping styles at genus level also support our previous results of relative abundance in Figures 2  and 3 as well as those of ternary plots in Figures S8-S10. Taxa that showed significant differential abundance are further analyzed in heat maps of the top 35 most abundant phyla and genera at the sample (Figures S15 and S16, respectively) and grouping style (Figures 7 and 8, respectively) levels.
Significant results from t-test between pairs of subgrouping style ABC of bulk soil (or group J), of subgrouping style DEF of rhizosphere soil (or group K) or across soil type (or grouping style GHI) are shown in Figures S17-S22 and referred to by arrows in Figures 7 and 8, respectively. The results in Figure 7 indicate significant differences in 11 phyla for grouping styles ABCDEF and JK. Of which, three phyla (e.g., Euryarchaeota, Verrucomicrobia and Nitrospirae) show significant differences in grouping style ABCDEF, five phyla (e.g., Cyanobacteria, Firmicutes, Chloroflexi, Actinobacteria and Proteobacteria) in grouping style JK, 12 and three phyla (e.g., Bacteroidetes, Acidobacteria and Gemmatimonadetes) in grouping style GHI. Among them, two phyla (e.g., Verrucomicrobia and Nitrospirae) have low relative abundance levels, thus, removed from further analysis. In terms of grouping style GHI, five phyla (e.g., Verrucomicrobia, Gemmatimonadetes, the archaeal phylum Thaumarchaeota, Acidobacteria and Euryarchaeota) show significant differences, of which, phyla Verrucomicrobia and Thaumarchaeota show low relative abundance levels, thus, removed from further analysis (Figure 7). Figure 7. Relative abundance for the top 35 taxa at the phylum level within grouping styles ABCDEF and GHI of microbiomes collected from surrounding bulk (grouping style ABC or group J) and rhizosphere (grouping style DEF or group K) soils of Nitrosalsola vermiculata after 0 (grouping style AD), 24 (grouping style BE) and 48 h (grouping style CF) of watering Samples data were also grouped based on time after watering regardless of soil type, e.g., 0 (group G), 24 (group H) and 48 h (group I). Arrows refer to phylum with significant differences within grouping styles. Arrows refer to the phyla that showed significant differences in terms of relative abundance by t-test. Black arrows refer to significant phyla at grouping style ABCDEF or at grouping style GHI, blue arrows refer to significant phyla at grouping style JK, while red arrows refer to significant phyla at the two grouping styles. Dotted arrows refer to phyla among the 10 most abundant ones (see Figure 2). Statistical analysis in terms of t-test bar plots is shown Figures S17-S22.
The results in Figure 8 indicate significant differences in 20 genera for grouping styles ABCDEF and JK. Of which, 12 genera (e.g., Massilia, Blastococcus, Rubellimicrobium, Planomicrobium, Flavisolibacter, Nocardioides, Ramlibacter, Devosia, Lysinibacillus, Domibacillus, Bryobacter and Streptomyces) show significant differences in grouping style ABCDEF, five genera (e.g., Bacillus, Ammoniphilus, Sphingomonas, Microvirga & unidentified genus of Cyanobacteria) in grouping style JK, and three genera (e.g., Pontibacter, Adhaeribacter and Arthrobacter) in grouping style GHI. Among them, all genera of grouping style ABCDEF show low relative abundance, thus, removed from further analysis. In terms of grouping style GHI, six genera (e.g., Streptomyces, Cellulomonus, Ramlibacter, Devosia, Flavisolibacter and Stenotrophobacter) show 13 significant differences. All these genera show low relative abundance, thus, removed from further analysis ( Figure 8). Therefore, a total of nine phyla (e.g., Cyanobacteria, Euryarchaeota, Chloroflexi, Actinobacteria, Proteobacteria, Bacteroidetes, Firmicutes, Acidobacteria and Gemmatimonadetes) and eight genera (e.g., Bacillus, Ammoniphilus, Sphingomonas, Microvirga, Pontibacter, Adhaeribacter, Arthrobacter and unidentified genus of Cyanobacteria) were analyzed further against the three grouping styles ABCDEF, GHI and JK. Figure 8. Relative abundance for the top 35 taxa at the genus level within grouping styles ABCDEF and GHI of microbiomes collected from surrounding bulk (grouping style ABC or group J) and rhizosphere (grouping style DEF or group K) soils of Nitrosalsola vermiculata after 0 (grouping style AD), 24 (grouping style BE) and 48 h (grouping style CF) of watering Samples data were also grouped based on time after watering regardless of soil type, e.g., 0 (group G), 24 (group H) and 48 h (group I). Arrows refer to phylum with significant differences within grouping styles. Arrows refer to the genera that showed significant differences in terms of relative abundance by t-test. Black arrows refer to significant genera at grouping style ABCDEF or at grouping style GHI, blue arrows refer to significant genera at grouping style JK, while red arrows refer to significant genera at the two grouping styles. Dotted arrows refer to phyla among the 10 most abundant ones (see Figure 2). Statistical analysis in terms of t-test bar plots is shown Figures S17-S22.
The results in grouping style ABCDEF indicate significant increase in relative abundance of phyla Acidobacteria and Gemmatimonadetes in rhizosphere soil at 24 and 48 h time points, while at 24 h only for phylum Bacteroidetes and 48 h only for the archaeal phylum Euryarchaeota ( Figure S17). At genus level, the results indicate significant increase in relative abundance of genera Pontibacter, Adhaeribacter and Arthrobacter in rhizosphere soil at 24 h time point, while genus Pontibacter also shows significant decrease at 48 h time point ( Figure S18). Relative abundance of phyla Acidobacteria, Gemmatimonadetes and archaeal Euryarchaeota in grouping style GHI shows significant increase at 24 h time point across soil types, while phylum Gemmatimonadetes also shows significant increase at 48 h time point ( Figure S19). However, none of the genera with high relative abundance show significant results in grouping style GHI ( Figure S20). In terms of grouping style JK, relative abundance of phylum Cyanobacteria in bulk soil is significantly higher than that in rhizosphere soil ( Figure S21). Opposite results are reached for the Gram-positive, e.g., Firmicutes and Actinobacteria, and Gram-negative, e.g., Proteobacteria Bacteroidetes, Acidobacteria, Chloroflexi and Gemmatimonadetes, phyla. Results in grouping style JK indicate significant higher relative abundance of unidentified genus of Cyanobacteria in bulk soil, while rhizosphere soil for genera Adhaeribacter and Pontibacter of phylum Bacteroidetes, Bacillus and Ammoniphilus of Firmicutes, Arthrobacter of Actinobacteria, Microvirga and Sphingomonas of Proteobacteria ( Figure S22).

Discussion
Differential abundance due to soil type In the present study, bulk soil refers to a native spot near the selected plants where the closest growing flora is at least 10 m 2 away. It is likely that high abundant bacteria in the plant rhizosphere interact mainly via symbiotic association, a phenomenon similar to legume-rhizobium symbiosis (Nilsson et al., 2005). Through this association, symbiotic bacteria receive main carbon source from the plant host, while, in return, bacteria provide nitrogen in several forms, ex., citrulline and glutamate (Lindblad et al., 1991;Bergman et al., 2007;Cuddy et al., 2012). Microbes can also promote plant growth by supplying substances such as IAA, amino acids, vitamins, solubilized phosphates and oxygen to root rhizosphere of the plant host (Liesack et al., 2000;Sinha and Häder, 2008;Prasanna et al., 2009;Liberton et al., 2013;Ranjan et al., 2016). Phylum Acidobacteria, for example, was reported to have strong symbiotic association with plant host (Kalam et al., 2020), where its enzymes metabolize, not only carbon, but also sulfur (Tank and Bryant, 2015;Eichorst et al., 2018).
Relative abundance of phylum Cyanobacteria or its descendent unidentified genus is significantly higher in bulk soil compared with that in rhizosphere soil (Figures 7 and 8). This indicates that this microbe might favors living independently in bulk soil although several reports indicated symbiotic associations with a wide range of plants (Liesack et al., 2000;Sinha and Häder, 2008;Liberton et al., 2013;Ranjan et al., 2016). In the present study, 22 OTUs of unidentified Cyanobacteria exist. Of which, OTU1 comprises the most abundant taxon across bulk and rhizosphere soils (Table S1). Interestingly, these OTUs are unidentified down the phylum level indicating that they represent new taxa of this microbe that live mainly in the bulk soil. Previous reports indicated that microbial abundance in rhizosphere soil can be as high as 10-1000 times that in bulk soil (Miransari, 2011;Zuo et al., 2021). The latter results meet those of the present study for rhizosphere soil of N. vermiculata, except for Cyanobacteria, where 19, out of the 22, OTUs are far higher in their relative abundance in bulk soil compared with that in rhizosphere soil (Tables 1 and S1).
In the present study, phyla Proteobacteria, Firmicutes, Acidobacteria, Actinobacteria, Bacteroidetes, Chloroflexi and Gemmatimonadetes show significantly higher microbial abundance in rhizosphere soil compared with that in bulk soil (Figures 2 and 7). Members of these phyla were reported to act as plant growthpromoting bacteria (Chukwuneme et al., 2020;Barnard et al., 2013;Yadav et al., 2018). Descendent genera of phyla Bacteroidetes (e.g., Adhaeribacter and Pontibacter), Firmicutes (e.g., Bacillus and Ammoniphilus), Actinobacteria (e.g., Arthrobacter), and Proteobacteria (e.g., Microvirga and Sphingomonas) in the present study also show significantly higher microbial abundance in rhizosphere soil compared with that in bulk soil (Figures 2 and 8). Genus Adhaeribacter (Bacteroidetes) was proposed to predominantly exist in rhizosphere of Broomcorn Millet (Na et al., 2019), while genus Pontibacter (Bacteroidetes) was reported to contribute to the hormetic responses of soil alkaline phosphatase (Fan et al., 2018). Bacillus (Firmicutes) is a highly tolerant genus that is potentially used in generating biofertilizers to diverse environmental stresses due to its ability to form spores resistant to desiccation (Barnard et al., 2013). Arthrobacter (Actinobacteria) was proven to be an important genus in the microflora that promotes plant growth especially under drought stress. Members of genus Arthrobacter are resistant to desiccation and can survive for a long time under starvation conditions (Chukwuneme et al., 2020). Abundance of genera Ammoniphilus (Firmicutes) and Microvirga (Proteobacteria) was high in Calotropis procera rhizosphere (Ramadan et al., 2021). Genus Microvirga also acts as plant growth promoting bacteria (PGPB) (Wang et al., 2017). Genus Sphingomonas (Proteobacteria) is more abundant in the present study in rhizosphere soil and is known to promote growth of Arabidopsis 15 thaliana, however, it can survive in bulk soil by taking energy from degrading organic pollutants .
Explanation for the high abundance and mutual dominance of the above-mentioned taxa in rhizosphere of N. vermiculata and other plant species is partially due to the intensity of specific root exudates the bacteria receive from the plant host (Lundberg et al., 2012;Dai et al., 2019;Zhang et al., 2019). Among other factors that might influence high abundance of microbes in rhizosphere soil of N. vermiculata, Lopes et al. indicated that range of soil pH implies shaping of soil microbiomes (Lopes et al., 2021). Chodak et al. demonstrated that microbes favoring high soil pH will grow better in rhizosphere soil, while those favoring low pH (≤ 5) grow better in bulk soil (Chodak et al., 2015). As soil pH in the western region of Saudi Arabia is estimated to be 6.1, this might partially explain the high abundance of members of these phyla in rhizosphere soil. We also speculate that natural co-existence of the highly abundant pairs of Gram-negative phyla, e.g., Acidobacteria and Proteobacteria, on one hand, and pairs of Firmicutes (Gram-positive) and Bacteroidetes (Gram-negative), on the other hand, might contribute to the high abundance of other microbes in rhizosphere of N. vermiculata. Wei and colleagues indicated that co-existing phyla Acidobacteria and Proteobacteria promote sustainability of other cohabiting phyla, like Actinobacteria, Chloroflexi and Gemmatimonadetes in the present study (Wei et al., 2017). Co-existence is a form of cell-cell communication, also known as bacterial coaggregation. The term coaggregation refers to specific recognition and adherence of different microbial species towards possible formation of biofilm. Soil surfaces are previously suggested to be a good setting for invisible multispecies biofilm formation to distribute antibiotic resistance genes among cohabiting microbes (Wang et al., 2015). Stevens et al. also indicated that bacterial coaggregation can promote the dynamic transmission of microorganisms among different environmental niches (Stevens et al., 2015), while high prevalence of synergism in multispecies biofilm formation of soil microbes was proven in another study (Ren et al., 2015). Highly abundant phyla also co-exist in the soil in a way to imply minimal competition for resources (Wei et al., 2017), where Firmicutes, for example, favors lipid nutrients in biogas reaction, whereas Bacteroidetes favors starch nutrients (Kampmann et al., 2012).
Influence of microbes on abundance of cohabiting ones can also be explained through the mechanism of quorum sensing (QS). The latter is a type of intercellular communication mechanism that is influenced by extracellular signaling molecules called autoinducers (AIs) (Rutherford and Bassler, 2012). These molecules mediate bacterial population density in the environment. These signaling molecules are encoded by a variety of bacterial genes that allow bacteria to act in synchrony. AIs also mediate bacterial sporulation, antibiotic production and biofilm formation. Gram-negative bacteria usually use acylated homoserine lactones as autoinducers, while Gram-positive bacteria use processed oligo-peptides to communicate (Ng and Bassler, 2009;Williams and Camara, 2009). It was suggested that bacterial AIs elicit specific responses from host plant (Miller and Bassler, 2001). When reached to a certain threshold, QS promotes signals within bacterial kingdom, on one hand, and between bacteria and the host plant, on the other hand. In return, plant responds to the signals by releasing exudates that promote growth performance of the cohabiting microbes.
Nonetheless, triggering QS in Cyanobacteria by specific lactones also occurs and results in the formation of cyanobacterial blooms and production of the toxin microcystin LR (Herrera and Echeverri, 2021). This action can take place during symbiotic association with plant rhizosphere. However, higher abundance of this microbe in the bulk soil still, in the present study, requires an explanation. One answer refers to the mutual abundance of the two most abundant genera in the present study, e.g., unidentified Cyanobacteria in the bulk soil versus genus Bacillus (Firmicutes) in the rhizosphere (Figure 9) and that of the phyla Cyanobacteria in bulk soil versus Firmicutes in rhizosphere soil ( Figure S23). We claim a possible role of genus Bacillus of Firmicutes in quenching QS of Cyanobacteria in the rhizospheres soil. Aligning with our claim, genus Bacillus was recently reported to act in blocking Cyanobacteria blooming, thus, inhibiting the release of the harmful toxin to the plant root (Wu et al., 2021). The latter authors indicated that Bacillus sp. strain S51107 exhibits a strong algicidal activity against the Cyanobacterium Microcystis aeruginosa. This activity is regulated by the signaling peptide NprX or pre-ComX encoded by comQ gene thatt triggers QS in Bacillus (Stanley and Lazazzera, 2005;Wu et al., 2021). To prove the possible action of pre-ComX in inducing QS in phylum Firmicutes, ex., genus Bacillus, we have tested the differential expression of the gene encoding pre-ComX, e.g., competence pheromone precursor or comQ gene, in metatranscriptomes collected from the two types of soil using real time PCR ( Figure S24). Sequence of this gene was retrieved from the NCBI for Bacillus subtilis subsp. subtilis str. 168 (NC_000964.3:3256008-3256907). The results of real time PCR aligned with those of differential abundance of genus Bacillus or phylum Firmicutes, where relative expression level of the gene increases in rhizosphere soil compared with that in bulk soil.
Differential abundance due to watering Watering of drought-stressed plants rapidly deviates natural progression of microbiome development, which can be eventually restored after soil is dried up. Xu et al. speculated that one week following watering is a maximum time to restore microbiome pattern in rhizosphere of drought-treated plant roots (Xu et al., 2018).
In the present study, two days after watering are enough to recover the natural growth pattern of some microbes at phylum and genus levels and to restore their previous or natural pattern under drought condition (Figures 7  and 8).
While, it seems that this short period is not enough for the majority of microflora to do both actions. This supports our claim that growth dynamics and genetic background of microbiota affect differential response to stress and, in some cases, can affect plant root architecture.
Generally speaking, Gram-positive (diderm) bacteria are known to be more drought tolerant than Gram-negative (monoderm) bacteria as the first has stronger cell wall and drought avoidance strategies (ex., spore formation) (Potts, 1994). Drought and rewetting stress were proven to significantly alter structure of soil metagenomes (Chodak et al., 2015). The latter authors indicated that abundance of Gram-positive bacterial phyla increases after the stress, whereas Gram-negative bacteria decreases. Nonetheless, the latter statement cannot be generalized as it seems that plant growth stage also affects relative abundance of rhizosphere microbes under drought stress where abundance of phylum Actinobacteria (Gram-positive) significantly increased at early stages in peanut, while that of phylum Cyanobacteria (Gram-negative) increased at flowering stage (Dai et al., 2019). Phyla Acidobacteria and Gemmatimonadetes of grouping style ABCDEF of the present study show higher relative abundance in rhizosphere soil at 24 and 48 h watering time points, while at 24 h only for phylum Bacteroidetes and 48 h only for archaeal phylum Euryarchaeota ( Figure S17). In addition, the bacterial phylum Acidobacteria and the archaeal phylum Euryarchaeota of grouping style GHI show higher relative abundance at 24 h watering time point regardless of soil type, whereas at 24 and 48 h for phylum Gemmatimonadetes ( Figure S19). At genus level, genera Pontibacter (Bacteroidetes), Adhaeribacter (Bacteroidetes) and Arthrobacter (Actinobacteria) in rhizosphere soil also show higher relative abundance at 24 h time point in rhizosphere soil, while genus Pontibacter only show further significant decrease at 48 h time point ( Figure S18). Phyla Acidobacteria, Gemmatimonadetes and Bacteroidetes are Gram-negative bacteria, while archaeal phylum Euryarchaeota can be either Gram-positive or Gram-negative based on the presence of pseudomurein in the cell wall (Iverson et al., 2012).
Highly abundant microbes whose abundance was not changed due to watering are the Gram-positive phylum Firmicutes and Gram-negative phyla Proteobacteria, Chloroflexi and Cyanobacteria and their descendent genera Bacillus, Ammoniphilus, Sphingomonas, Microvirga and unidentified genus of Cyanobacteria. We speculate taxa that do not respond to the sudden presence of water and keep their growth rate regardless of the change in water supply are drought tolerant. Firmicutes and Proteobacteria are considered in several reports to have potential significance in assisting plants to withstand drought stress (Jang et al., 2020;Ma et al., 2020). Jang et al. indicated that relative abundance of phylum Proteobacteria in bulk soil was lower under drought condition than that under watered conditions and vice versa in the rhizosphere soil of rice (Jang et al., 2020). Bacterial responses to drought and rewetting stress were shown to be based on genetic makeup of the bacteria, host-microbe interaction as well as nutrient supply in the surrounding environment, e.g., contents of soil N and organic C, heavy metal, etc. (Chodak et al., 2015). Relative abundance of Actinobacteria under drought stress was shown in one study to be higher in the rhizosphere than that in the surrounding bulk soil (Naylor et al., 2017) as this phylum was proven to survive under unfavorable environmental conditions (Passari et al., 2015). The phylum is also known to produce bioactive secondary metabolites acting as plant growth promoting (PGP) and nitrogenous compounds in non-legumes and to make P solubilization via production of various organic acids including citric, gluconic, lactic, malic, oxalic, propionic and succinic acids (Sathya et al., 2017;Yadav et al., 2018). Nonetheless, Breitkreuz and colleagues indicated that relative abundance of Proteobacteria decreased under drought stress in wheat rhizosphere, while increased for Chloroflexi and Firmicutes (Breitkreuz et al., 2021). Two other studies on rice and cotton rhizosphere microbiomes revealed that phylum Chloroflexi dominates under drought conditions (Santos-Medellín et al., 2017;Ullah et al., 2019).
Thus, the four phyla Firmicutes, Proteobacteria, Actinobacteria and Chloroflexi are considered, in several reports, to have potential significance in assisting plants to withstand drought stress (Ma et al., 2020). These results are consistent with those of many other reports on many other plants (Chodak et al., 2015;Bu et al., 2018;Gumiere et al., 2019).
A bacterial strain of Bacillus sp. was proven to produce and secrete spermidine (SPD) in Arabidopsis thaliana (Zhou et al., 2016). SPD is a type of polyamine (PA) that plays a role in plant growth under drought stress. The authors indicated that inoculation with the bacteria increases plant biomass, improved root system architecture and plants have displayed stronger ability to tolerate drought stress than non-inoculated plants. This ability is due to the higher levels of abscisic acid (ABA) under drought stress. In addition, stress tolerant Bacillus was reported to induce physiological response in the plant root vicinity in order to help the plant alleviating adverse effects of drought stress (Bokhari et al., 2019). Genus Sphingomonas is known for its ability to tolerate drought stress and promote tolerance to the neighboring plant roots (Luo et al., 2020). Earlier report of the latter research group indicated that a species of genus Sphingomonas was able to promote the growth of A. thaliana by driving roots plasticity during development, on one hand, and inducing changes in bacterial community in the plant root rhizosphere, on the other hand. The same study indicated high abundance of genus Sphingomonas under drought stress, suggesting the important role of the microbe in conferring drought 18 stress tolerance to the plant. Up to our knowledge, there are no previous reports in the literature referring to the role of genera Ammoniphilus and Microvirga under drought stress.

Conclusions
In the present study, differential microbiota signatures were detected in the rhizosphere of the wild plant A. vemiculata and surrounding bulk soil. We justified differential abundance of microbes in the two soil types by the involvement of inter-microbial as well as plant-microbe interactions, on one hand, and microbial coexistence as well as quorum sensing that refers to the dynamics of microbial intensity in different environmental niches, on the other hand. Upon watering, some soil microbes respond positively, and the 48 h was enough for some of these microbes to restore their natural growth pattern under drought stress. Nonetheless, nonresponsive microbiota was considered as drought tolerant. In conclusion, the present study warrants further research to dissect factors involved in shaping microbiome signature in different soil types and factors influencing differential stress responses of microbiota. The study also recommends harnessing the rhizosphere microbiomes to increase the resilience of crop production to drought stress conditions in the future.