The Z. qatarensis plants (shoot and root), rhizosphere and bulk soil were collected from four major population locations (T1, T2, T3, and T5) in the extreme desert of Empty Quarter in Oman (Supplementary Table 1 and Figure 1). The Z. qatarensis plant populations are distributed in a conserved and consistent ecological pattern across the four distinctive locations in desert areas. Each population was approximately 100 km apart and samples within a population were collected in triplicate. Each replica from the individual population was representative of ten plants/soils that were pooled for DNA extractions and chemical analysis. The rhizosphere soils adjacent to the root surface (10 to 15cm deep) were collected (Supplementary Figure 1). These were mature plants with a reasonably defined rooting system and the sandy soils were removed with the help of a sterile blade. There were no specific signs of rhizosheath with the root, so the soil attached to root parts was categorized as rhizosphere – following the classification of Pang et al. (2021). The bulk soil samples were collected from a depth of 10 to 15cm with no signs of Z. qatarensis presence. The root parts were carefully collected by removing the sand and particulate matters. For endosphere, the root and shoot samples were washed with sterile distilled water and sodium hypochlorite to remove epiphytic microbes following the method of McPherson et al. (2018). Briefly, the pre-sterilized scalpel (with 70% ethanol) was used to prepare individual tissues from about 10 to 15 roots and shoot parts ranging from 4 to 6 cm in length. The excised tissues were placed in autoclaved phosphate buffer (NaH₂PO₄ − 6.33 g/L, Na₂HPO₄ − 8.5 g/L with pH = 6.5: McPherson et al. (2018). The sterilized samples were stored at − 20°C for DNA extractions. All the samples were stored at 4°C for soil chemical analysis. The sampling area’s climate is dry with annual rainfall lower than 30 mm and summer temperature can reach up to +48°C with a relative humidity of 20–30% (Supplementary Table 2). The samples were collected during the dry summer season (May–June 2020).

After pooling the soil/plant samples within each replicate, 10.0 g mixtures were subjected to total DNA extraction using the MoBio Power Soil DNA Extraction Kit. The DNA was quantified using Qubit 4.0 and high sensitivity kit (Invitrogen, United States). Libraries of each DNA sample were generated by amplifying the internal transcribed spacer (ITS region) and 16S rRNA gene (V4 region) using Nextflex PCR I Primer Mix (Perkin Elmer, United States) for fungal and bacterial communities, respectively (Supplementary Figure 2). Illumina amplicon primers for 16S (Forward 5′–TCGTCGGCAGCGTCAGATGTGTATAAGAGACAGCCTACG GGNGGCWGCAG–3′; Reverse 5′– GTCTCGTGGGCTCG GAGATGTGTATAAGAGACAGGACTACHVGGGTATC TAA TCC –3′) and ITS (Forward 5′–AATGATACGGCGACCAC CGAGATCTACACGG CTTGGTCATTTAGAGGAAGTAA–3′; Reverse 5′–CAAGCAGAAGACGGCATACGAGA TCGGCTG CGTTCTTCATCGATGC–3′) were used. For the 16S rRNA gene, peptide nucleic acid (PNA) clamps were used to reduce mitochondrial and chloroplast contamination.

A paired-end sequencing approach with read lengths of 250 bp was conducted on an Illumina MiSeq instrument (Illumina Inc., San Diego, CA, United States) operating with v2 chemistry (User Guide Part # 15027617 Rev. L). All quality reads related to the study are available at NCBI under BioProject PRJNA771947 (SRP341951) and PRJNA767523 (SRP339516) for bacteria and fungi, respectively.

The sequencing reads were analyzed with QIIME 2.0 (Bolyen et al., 2019). First, reads from ITS and 16S rRNA amplicons were separated into different files. Then, the average quality of forward and reverse reads was observed in each dataset. Only the forward reads were used for the following analyses due to the low quality of the reverse reads in both datasets. We used the DADA2 algorithm for denoising and generating the amplicon sequence variants (ASV) table for each dataset (Callahan et al., 2016). The reads clustered in the same ASV have nucleotide sequences that are 100% identical (Callahan et al., 2017). In the denoising, sequences were filtered by overall quality, trimmed in low-quality regions, and chimeric sequences were removed. The 16S rRNA gene reads were trained on the SILVA database for taxonomic classification (Quast et al., 2012), while the UNITE database was used to classify the ITS sequences (Nilsson et al., 2019). Sequences classified as mitochondria and chloroplast were removed from the 16S rRNA gene ASV table. The 16S rRNA gene and ITS ASV tables were rarefied to 6000 and 100 reads for diversity analyses containing all samples from each dataset, respectively. The ASV tables of each dataset were then split for each sampling compartment for a more detailed analysis. For the bacterial diversity analyses, the bulk soil, rhizosphere, root endosphere, and shoot endosphere ASV tables were rarefied to 30,000, 37,000, 7,500, and 6,200 reads, respectively. For the fungal diversity analyses, the bulk soil, rhizosphere, root endosphere, and shoot endosphere ASV tables were rarefied to 120, 100, 1,300, and 1,200 reads, respectively. For beta-diversity analyses, the Bray-Curtis distance matrix was generated for each dataset and exported to statistical software. For alpha-diversity analyses, the Shannon index and Observed ASVs were calculated in each dataset to infer species diversity and richness, respectively.

The Bray-Curtis distance matrices were exported to R for multivariate statistical analyses. Non-metric multidimensional scaling (NMDS) was performed with the metaMDS function to show the ordination of samples according to the major factors (plant compartment) affecting the microbial communities, while constrained analysis of principal coordinates (CAP) was performed with the capscale function to show the ordination of samples according to more specific factors (sampling locations and soil type) using the “vegan” package v. 2.5-6 (Dixon, 2003). Both NMDS and CAP results were visualized with the ggplot function in the “ggplot” package v. 3.3.0 (Wickham, 2011). Permutative multivariate analysis of variance (PERMANOVA, 999 permutations) was used to test for significant effects of the factors studied (plant compartment, sampling location and soil type) and their interaction on bacterial and fungal beta-diversity using the “adonis function.” Analysis of similarity (ANOSIM) was also used as a second method to test the effects of the factors studied on the fungal and bacterial communities using in R. Differences in species diversity (Shannon index) and richness (Observed ASVs) for the same factors were assessed using the Kruskal-Wallis test in QIIME 2.0 (Bolyen et al., 2019). Changes in the relative abundance of microbial genera between soil types were evaluated with the Kruskal-Wallis test and Bonferroni p-value correction using STAMP software (Parks et al., 2014).

The four populations were sampled from soils characterized as sandy (T1/T3) and gravel (T2/T5), where the pH ranges from 6.7 to 7.2 and electrical conductivity from 27 to 39 dS m^–1. The soils with more gravel structure have higher nitrates than sandy ones (Supplementary Table 2). A detailed soil chemical analysis was performed (Khan et al., 2020). Furthermore, essential nutrients such as silicon (Si), magnesium (Mg) and calcium (Ca) were quantified from plants and soil from the four populations as described previously (Bilal et al., 2018) using inductively coupled plasma mass spectrometry (ICP-MS; Optima 7900DV, Perkin-Elmer, United States). All the measurements were carried out in triplicate.

At least three replicate samples were analyzed during this study. The data for the enzyme study is presented as the mean ± standard error (SEM). The significant differences were determined using one-way analysis of variance (ANOVA). The differences were considered significant at P < 0.05 and were calculated by GraphPad Prism Version 6.01 (GraphPad Software, San Diego, CA, United States). Duncan’s multiple range test at P < 0.05 (SAS 9.1, Cary, NC, United States) was used to compare the mean values.

A total of 58,648 ITS and 2,715,959 16S rRNA sequence reads passed all quality filters and were used for the following analyses. We analyzed the differences in bacterial/archaeal and fungal community structure between four locations where Z. qatarensis was naturally found (Figure 1A). Location T1 and T3 have sandy soils, while locations T2 and T5 have gravel soils (Figure 1A). The microbiomes of two soil (bulk soil and rhizosphere) and two plant endosphere (roots and shoots) compartments were investigated (Figure 1B and Supplementary Figure 3). The bacterial/archaeal community was dominated by the phylum Proteobacteria in the root and shoot endosphere compartments and Acidobacteria in the rhizosphere and bulk soil compartments (Figure 2A). Proteobacteria had a higher relative abundance in sandy than gravel soils in the bulk soil and shoot endosphere. In the rhizosphere and root endosphere compartments, the phylum Firmicutes was relatively more abundant in sandy than in gravel soils (Figure 2A). On the other hand, in the shoot endosphere Firmicutes were proportionally more abundant in gravel than in sandy soils (Figure 2A). PERMANOVA results indicated that sampling compartment was the main factor significantly affecting (p < 0.001; R² = 0.34) the bacterial/archaeal communities, followed by location (p < 0.001; R² = 0.08). The interaction effect between compartment and location was also significant (p = 0.031; R² = 0.12) on the bacterial communities. The bulk soil and rhizosphere samples did not show a clear separation in the NMDS, but clustered separately from the root and shoot endosphere samples (Figure 2B). ANOSIM results also showed a significant effect of compartments (p < 0.001; R = 0.53) and locations (p = 0.022; R = 0.34) on the bacterial communities. When analyzing the samples from all compartments together, soil type was not a significant factor in shaping the bacterial community structure according to PERMANOVA and ANOSIM.

Bacterial alpha-diversity was also significantly different between compartments, with greater species richness and diversity in the rhizosphere, followed by bulk soil, root endosphere, and shoot endosphere (Figures 2C,D). The bacterial community distribution and abundances across different soil types and samples have also been shown through heat and networking maps (Supplementary Figures 4, 5). The analysis of each compartment separately indicated that the bacterial community structure was significantly different between the gravel and sandy soils in the rhizosphere (p = 0.025; R² = 0.07) and bulk soil (p = 0.013; R² = 0.09), but locations had a larger effect than soil type for both rhizosphere (p < 0.001; R² = 0.29) and bulk soil (p = 0.02; R² = 0.17) (Figure 3). Location T2 had the most different bacterial community composition among the four locations (Figure 2). The root and shoot endosphere bacterial community structure was not affected by soil type and location. The alpha diversity was not affected by soil type and locations in any of the four compartments.

The Ascomycota phylum dominated the fungal community in all samples from the four compartments (Figure 4A). This dominance was mainly represented by the ascomycete species Acidea extrema. The root and shoot endosphere compartments had proportionally more unclassified sequences than the bulk soil and rhizosphere (Figure 4A and Supplementary Figure 6). As observed for the bacterial community, the fungal community structure was more affected by plant compartment (p < 0.001; R² = 0.21), followed by sampling location (p = 0.014; R² = 0.08) according to PERMANOVA. The interaction between these two factors also significantly affected the fungal communities (p = 0.032; R² = 0.16). However, ANOSIM indicated that only location (p = 0.028; R = 0.05), but not compartment (p = 0.119; R = 0.05) significantly affected the fungal community structure. The separation between samples from the four compartments was unclear in the NMDS ordination (Figure 4B). Soil type did not affect the fungal beta-diversity when analyzing the samples from the four compartments together.

As seen for bacteria, the fungal alpha-diversity was also affected by plant compartments, with greater species richness and diversity in the rhizosphere and bulk soil than in the root and shoot endosphere (Figures 4C,D). It is noteworthy the low fungal alpha diversity in all compartments. Rarefaction curves showed that the sequencing depth and sample rarefaction (100 reads) used were sufficient to reflect the fungal diversity in these samples (Supplementary Figure 7). In addition, the species richness was also affected by different locations, with greater values in T3 and T5 than in T1 and T2 (Supplementary Figure 5). When separating the analysis by each compartment, soil type still did not affect the fungal beta-diversity in any compartment. Soil type also did not affect the fungal alpha-diversity in any compartment. However, location significantly affected the fungal community structure of the rhizosphere (p < 0.001; R² = 0.37) and bulk soil (p = 0.049; R² = 0.23) compartments (Figure 5). Like what was observed for bacteria, the fungal communities of location T2 were the most different among the four locations, mainly in the rhizosphere (Figure 5). Moreover, the fungal species diversity decreased in T2 compared to T3 and T5 in the rhizosphere compartment (Supplementary Figure 5).

After analyzing the microbiome differences between locations and soil types, we aimed to identify the core microbiome of Z. qatarensis regardless of soil type and geographical location. We considered the core microbiome to be the microbes shared by >90% of all plant samples along with the four locations. Results indicated that no specific bacterial ASV was in the core microbiome of the bulk soil, rhizosphere, and root endosphere compartments, but one ASV classified in the Alcaligenes genus was present in the shoot endosphere of all plant samples. Similarly, no fungal ASV was present in the bulk soil and rhizosphere core microbiomes. However, two ASVs were found in the root endosphere core microbiome and three ASVs in the shoot endosphere core microbiome. One of these ASVs was in the core microbiome of both root and shoot endosphere and was classified in the species Acidea extrema, while the others were from unknown fungal genera. Next, we analyzed the core microbiome at the genus level. A total of 27 bacterial genera were in the bulk soil core microbiome, including Bryobacter, Chthoniobacter, Gemmata, Nitrospira, Pirellula, Rubrobacter and other 21 unknown/undescribed genera. On the other hand, no bacterial genus was in the rhizosphere core microbiome, while only one genus was in the root and shoot endosphere core microbiomes: Alcaligenes. Only the fungal genus Acidea (species A. extrema) was in the core microbiomes of all compartments. These results indicate that the different locations and soil types significantly affect the microbial beta-diversity and the core microbiome – shrinking it to none or just a few microbes.

We further investigated the changes in the microbial communities between the soil types by analyzing shifts in the relative abundance of specific bacterial and fungal genera. The bulk soil compartment had the highest number of bacterial genera (30) with significantly different relative abundances between the gravel and sandy soils, from which 18 were enriched in the sandy soils and 12 were enriched in the gravel soils (Figure 6). The bacterial genus with the largest relative abundance in sandy soils compared to gravel soils was Adhaeribacter. The genus with the largest relative abundance in gravel compared to sandy soils was Rubrobacter (Figure 6). The rhizosphere was the second compartment with more changes in the relative abundance of specific bacterial genera (17) between the soil types, where 13 genera were enriched in the sandy soils and four genera were enriched in the gravel soils (Figure 6). The bacterial genus with the most extensive proportional enrichment in the sandy soils was Bacillus. Blastocatella was the bacterial genus with the most significant relative abundance in gravel soils compared to sandy soils (Figure 6). The root endosphere showed only three changes in the relative abundance of bacterial genera between the soil types, including the enrichment of Streptomyces and Cupriavidus in the sandy soils and the putative genus WD101_soil_group in the gravel soils (Figure 6).

The genera affected by soil types were generally not the same in the rhizosphere and bulk soil, except Noviherbaspirillum, Marine_Group_II, and Herpetosiphon, which were enriched in the sandy soils in both compartments (Figure 6). Furthermore, Streptomyces was enriched in sandy soils in both rhizosphere and root endosphere (Figure 6). No bacterial genera had different relative abundances between soil types in the shoot endosphere compartment. No fungal genera changed in relative abundance between soil types in any of the four compartments. Many bacterial genera also had differences in relative abundance between locations within each soil type, but there were more changes in the rhizosphere than in the bulk soil (Supplementary Figures 4, 5). The root and shoot endosphere did not show any differences in genera between locations, and no fungal genera differences in relative abundance between locations in any compartment.

Three major nutrients, i.e., silicon (Si), magnesium (Mg), and calcium (Ca), were quantified in the soil, rhizosphere, and shoots across four populations of Z. qatarensis. In addition to major nutrients, we also assessed the soil physical and morphological properties shown in Supplementary Figure 8. The results showed that the rhizosphere had significantly (p < 0.01) higher Si content than bulk soil. Among the populations growing in each soil type, T3 had a significantly higher Si content than T1 in sandy soils, whereas T2 had a significantly higher Si content than T5 in gravel soils (Figure 7). The rhizosphere of T1, T3, and T2 had 23, 21.5, and 15.8% more Si than bulk soil. However, the rhizosphere of T5 had a significantly lower Si content compared to bulk soil. In the case of root to shoot parts, T1 and T3 (sandy soils) had significantly higher (p < 0.05; 8.7 and 41.6%, respectively) Si content in roots compared to T2 and T5 (gravel soils). On the other hand, the T2 and T5 had significantly higher (p < 0.05; 13.2 and 11.8%, respectively) Si in shoots compared to roots (Figure 7).

In the case of Mg, T2 (gravel texture) had a significantly higher (p < 0.003; 44.7%) content in bulk soil, while T5 had a higher content in rhizosphere. In contrast, T1 and T3 (sandy) had significantly higher (p < 0.001; 4.2 and 21.8% respectively) Mg content in the rhizosphere than in bulk soil. Overall, Mg amounts were more elevated in gravel than in sandy soils across the different populations. In plants, the Mg content was 7- to 8-fold higher (p < 0.001) in shoots than in roots. Among populations, shoots of T1 and T5 had a significantly higher Mg content than T3 and T2. Only T5 had a higher Mg content than T1, T2, and T3 (Figure 7).

In the case of Ca, it was significantly higher in the rhizosphere of T1 (8.7%) and T3 (7.2%) compared to the sandy bulk soils (T2 and T5) of Z. qatarensis. Conversely, in gravel soil, the T2 population of Z. qatarensis had a significantly higher (14.2%) Ca content in bulk soil compared to rhizosphere. T5 had lower Ca than the other populations of Z. qatarensis growing in different locations. In the case of the plants’ organs, shoots of Z. qatarensis had in general, a significantly higher Ca content than roots. A considerably higher amount of Ca was observed in T2 > T1 > T3 > T5 in shoots. This increase was four to fivefold higher than in roots (Figure 7).