We geochemically defined seventeen space-for-time soil transects, representing a wetness gradient from the shores of six lakes during the 2016 and 2017 Antarctic field seasons. The lakes span the length of the Wright (Lake Brownworth and Lake Vanda) and Taylor (Spaulding Pond, Lake Fryxell, Lake Hoare, and Lake Bonney) valleys (Figure 1 and Table 1). Lake Brownworth and Spaulding Pond are open-basin lakes, meaning that they have an outlet channel through which water flows out. The presence of the outlet regulates the water level in these lakes (Levy et al., 2018). Transects defined from these lakes are referred to in this study as static (as having more static water levels). Lakes Vanda, Fryxell, Bonney, and Hoare are closed-basin lakes, which do not have an outlet channel for water to flow out of the lake. The water levels in these lakes are dynamic on multiple time scales but, over the last 100 years, are thought to have been intermittently rising (Castendyk et al., 2016; Levy et al., 2018; Ramoneda Massague et al., 2021). Transects defined from these lakes are referred to in this study as dynamic (as having more dynamic water levels). For each lake, except for Lake Hoare, three space-for-time transects were defined along the wetness gradient. At Lake Hoare, only two transects were sampled due to snowfall while sampling.

In these soils, pH, electrical conductivity (EC) and water activity (WA) have been described in several studies as being the key drivers of community assembly in the MDV (Barrett et al., 2009; Van Horn et al., 2014; Niederberger et al., 2015, 2019; Feeser et al., 2018; Bottos et al., 2020; George et al., 2021). These metrics were measured in the field to identify in situ geochemical gradients across transects from all water bodies within which three geochemical zones were identified (Table 1 and Supplementary Figures 1, 2). This process involved performing measurements with an average spacing of 50 cm from the existing water level within a transition zone to observe the natural variation of these geochemical variables along the wetness gradient (Supplementary Figures 1, 2). The first zone represented a water-saturated wet zone characterized by comparatively higher water activity (average of 1) and the lowest electrical conductivity (138 ± 81 μS) (sampling point 1). With increasing distance from the lake shore, water content remained high, but conductivity increased several times (938 ± 1,591 μS) (sampling points 2 and 3). This was used to define the second zone, representative of a transition zone. The third zone represented a dry zone, characterized by the lowest measured water activity and lower electrical conductivity (341 ± 443 μS) (sampling points 4 and 5).

Briefly, electrical conductivity and pH were determined from a 1:5 slurry of soil to Milli-Q water (Thermo Scientific Orion meter, United States) (Lee et al., 2012). Water activity was measured using a PawKit meter (AquaLab, NZ), following the manufacturer’s instructions. Infield measurements of electrical conductivity and pH were confirmed under laboratory conditions using Thermo Scientific Orion meter and the same slurry technique, but with the addition of 0.5 mL of 0.01M CaCl₂ per sample (Minasny et al., 2011; Lee et al., 2012). Soil moisture content was determined according to Lee et al. (2012). The elevation of each sampling point along the transects was surveyed in relation to the lakeshore (sampling point 1) using a laser leveling system (Topcon Laser Systems, United States).

For microbial community analysis, we aseptically collected surface soil samples (top 2 cm) with a sterile spatula from the five geochemically predefined sampling points. Spatulas were washed and then sterilized with ethanol wipes between sampling points and sampling sites. In total, across the six lakes, seventeen space-for-time transects were sampled, resulting in 85 soil samples. Individual samples were thoroughly mixed in a sterile Whirl-Pak^® bag and then distributed into sterile 50 mL Falcon tubes. Samples were initially stored in the field on ice for 24 h before being transferred to dry ice and transported to Scott Base (the New Zealand Antarctic Research Station), where they were maintained a –60°C until DNA extraction. The remaining samples in each Whirl-Pak^® bag were used to confirm moisture content, electrical conductivity, and pH analysis under laboratory conditions. These samples were kept cold at 4°C until processed.

Total DNA was extracted from approximately 1 g of soil using a modified version of the CTAB (cetyl-trimethyl-ammonium-bromide) bead-beating method (Coyne et al., 2001). For each batch of DNA extractions, a negative control was included to ensure the identification of possible contamination. DNA concentration was determined using the Qubit dsDNA HS Assay Kit (Thermo Fisher Scientific, United States). The bacterial and archaeal microbial community was targeted through the V4 region of the 16S rRNA gene was amplified in triplicate using the fusion-primer set 515F/806R (Parada et al., 2016) and sequenced using Ion PGM chemistry. Briefly, 20 μL PCR reactions each contained: 1 ng of total DNA, dNTPs (240 μM), MgCl₂ (6 mM), bovine serum albumin (0.24 μM), forward and reverse fusion primers (0.2 μM), 1 U of Platinum Taq polymerase (Invitrogen Inc., Carlsbad, California) and PCR buffer (1.2×). The following PCR conditions were used: 94°C for 3 min, followed by 30 cycles at 94°C for 45 s, 50°C for 1 min, 72°C for 1.5 min, and a final extension step at 72°C for 10 min. For each PCR run, negative amplification was confirmed for each DNA extraction batch control. Triplicate PCR amplicons were pooled, and the expected amplicon size was confirmed via electrophoresis with a 1% agarose TAE gel. PCR products were cleaned and each sample concentration was normalized with SequalPrep™ (Thermo Fisher Scientific, United States). Normalized samples were pooled at an equimolar concentration into a single library for sequencing. Amplicon sequencing was performed using the Ion PGM™ System for Next Generation Sequencing (Thermo Fisher Scientific, United States) at the Waikato DNA Sequencing Facility (University of Waikato, New Zealand).

Raw sequences were filtered with Ion PGM™ software to remove low-quality and polyclonal reads. The remaining sequences were processed using a combination of Mothur (v.1.40.5) and USEARCH 10 (v10.0.24) software (Schloss et al., 2009; Edgar, 2010, 2013). Forward and reverse primers were identified within the sequences and trimmed using the python script fastq_strip_barcode_relabel.py supplied by UPARSE (v10.0.240). Sequences without forward and reverse primers were discarded. The remaining sequences were trimmed based on the length and the number of homopolymers sourced by Mothur script (Schloss et al., 2009). All reads with expected error rates higher than 2.5 were discarded using USEARCH (Edgar, 2010), and all reads were truncated to 350 bp. Through dereplication, unique sequences were identified and their abundances quantified. Sequences were sorted, singletons removed, and the remaining sequences were clustered into representative OTUs using the UPARSE-OTU algorithm combined with the GOLD database to detect and remove chimeras (Bernal et al., 2001; Edgar, 2013). Reads were finally clustered into operational taxonomic units (OTUs) using UCLUST with a similarity threshold of 97%. Sequences that did not map to any OTU were discarded. Taxonomy was inferred using SINA (v1.2.11) and the SILVA SSU database (v 138) (Pruesse et al., 2012). A 0.005% cutoff was applied to the raw OTU counts across the dataset to remove poorly represented OTUs (Bokulich et al., 2013). Filtered OTUs were subsequently removed from the sequence Fasta file. Sequences were aligned using a Multiple Alignment using Fast Fourier Transform with default settings (MAFFT v7.429-gimkl-2020a) and a phylogenetic tree was generated using FastTree (v2.1.11).

All statistical and visualization analyses were computed in R (v3.5.2, R Core Team, 2000) using the following packages: phyloseq (v1.26.1) (McMurdie and Holmes, 2013), vegan (v2.5) (Oksanen et al., 2012), picante (v1.8) (Kembel et al., 2010), ggplot2 (v3.2.1) (Wickham, 2016), ggpubr (v0.4.0), compositions (van den Boogaart and Tolosana-Delgado, 2008), and randomForest (v4.6) (Liaw and Wiener, 2002). Alpha diversity was calculated on a non-rarefied dataset using richness and phylogenetic diversity (PD) indexes. Differences in library sizes were tested on both rarefied and non-rarefied data with no significant impact on the data interpretation. Diversity differences between transect zones and lakes were tested using one-way ANOVA. Normality assumptions were tested using a Shapiro-Wilk test, and the assumption of homogeneity of variances was tested using a Levene’s test. For the relationship between the microbial diversity and the log transformed moisture content and electrical conductivity, correlation analyses (Pearson) were applied. For beta diversity analysis we transformed the data using a total sum normalization. Beta diversity was calculated using weighted UniFrac phylogenetic pairwise distances (Lozupone and Knight, 2015) and visualized in a principal coordinates analysis (PCoA). A PERMANOVA analysis, using Adonis function in vegan package (Oksanen et al., 2012), tested the dissimilarities among communities from different groups (lakes, and zones within the wetness gradients). The variation between samples from the three wetness zones described along the transects was tested using a permutational multivariate analyses of dispersion (betadisp) from vegan package (Anderson et al., 2006). A Random Forest model (Liaw and Wiener, 2002) was used to detect the most reliable and relevant top 20 OTUs to predict the different zones along the wetness gradients. Raw counts were transformed using a centered log-ratio (clr) before classification and regression models. OTUs were selected based on the MeanDecreaseAccuracy values, which measures the extent to which a variable (OTU) improves the accuracy of the forest in predicting the classification. Higher values indicate that the OTU improves the prediction. To evaluate the phylogenetic community assembly among species within each sample we calculated the mean nearest taxon distance (MNTD) and the nearest taxon index (NTI) using “taxa.labels” null model, with 999 iterations, using the function “ses.mntd” from the R package picante. The NTI was quantified as the number of standard deviations that the observed MNTD was from the mean of the MNTD null distribution, multiplying by –1. For a single community, observed NTI values > + 2 or < −2 indicate phylogenetic clustering of species or phylogenetic overdispersion, respectively. NTI values between −2 and 2 usually indicates the influence of a stochastic assembly (Stegen et al., 2012). To compare community assembly processes along the moisture gradients, we calculated the Beta Nearest Taxon Index (βNTI). Following Stegen et al. (2013), the βNTI is the number of standard deviations that the observed beta mean nearest taxon distance (βMNTD) is from the mean of the null distribution. It indicates how much the observed difference between a pair of communities differs from a null distribution. Similarly to NTI, values > + 2 or < −2 indicate phylogenetic clustering of species or phylogenetic overdispersion, respectively, while values between −2 and 2 usually indicates the influence of a stochastic assembly (Stegen et al., 2013). β-NTI pairwise comparisons were plotted against increasing differences in moisture content. A Euclidean distance matrix was calculated using pH, electrical conductivity, water activity, moisture content, and elevation log-transformed and normalized data. ANOSIM was performed on the resemblance matrix to test the significance of the dissimilarities between the predefined zones of the wetness gradients. Sampling locations were plotted using Quantarctica (v3.1) (Matsuoka et al., 2021).

Seventeen geochemically defined space-for-time transects were parameterized using measurements of elevation, soil moisture content, electrical conductivity, pH, and water activity from both dynamic (n = 11 transects; Lake Vanda, Lake Bonney, Lake Hoare, and Lake Fryxell) and static (n = 6 transects; Lake Brownworth and Spaulding Pond) lakes across the Wright and Taylor Valleys (Figure 1, Table 1, and Supplementary Figure 2). Across all transects, elevation increased from wet (sampling point 1) to dry (sampling points 4 and 5) soils (from 0 to 135 ± 19 cm), and soil moisture content decreased in the same direction (from 17.25 ± 3.07 to 0.54 ± 0.54%) (Figure 2 and Table 1). Electrical conductivity was lowest in wet and dry soil samples (138 ± 81 μS and 341 ± 443 μS, respectively), achieving the highest values and variation within the transition zone (938 ± 1,591 μS) (Figure 2 and Table 1). Water activity decreased from wet to dry soils in Lake Vanda, Lake Bonney, Lake Brownworth, and Lake Hoare and remained stable across the wetness gradient in Spaulding Pond and Lake Fryxell (Figure 2 and Table 1). The variability in soil pH across the wetness gradients was lake-specific (Figure 2 and Table 1). The significance of the distinction between the three wetness zones (wet, transition, and dry) was confirmed with non-parametric testing (ANOSIM, R = 0.69, p-value < 0.01, 999 permutations) (Supplementary Figure 3).

After filtering out low-quality and short sequence reads, we obtained a total of 2,084,605 sequence reads and 4,098 OTUs at 97% sequence similarity (Supplementary Table 1). We then filtered the low abundant OTUs using a 0.005% relative abundance cutoff across the entire OTU dataset and removed 72% of the initial OTUs. The remaining 28% (1,133 OTUs) comprised 97% of the initial reads. A significant correlation exists between the ordinations of the initial and filtered datasets (m12 squared: < 0.001; Procrustes Correlation: 0.99; p-value < 0.01, 999 permutations). We removed thirty-two OTUs classified as Chloroplast using the SILVA database, leaving a dataset containing 1,101 OTUs for downstream analysis. Sixteen percent of the remaining OTUs were unclassified at the phylum level. The Procrustes Correlation between the ordination matrices with and without the unclassified phyla was significant (m12 squared: < 0.001; Procrustes Correlation: 0.99; p-value < 0.01, 999 permutations). As a result, we retained the unclassified OTUs in the analysis.

We observed significant differences in species phylogenetic diversity (PD) within the dry, transition, and wet transect zones from each lake (Supplementary Figure 4; p-value < 0.01). Spaulding Pond consistently presented the most diverse community along the wetness gradients and Lake Vanda the least diverse community (Supplementary Figure 4). Except for Lake Vanda and Lake Bonney, PD correlated positively with soil moisture content (p-value < 0.05; Figure 3A). Soil electrical conductivity correlated negatively with PD (p-value < 0.05; Figure 3B) only in Lake Bonney transects, likely driven by the comparatively high electrical conductivity (5,300 μS) measurements.

The soil microbial community was represented by 24 phyla, with 98% of reads classified within the top 10 phyla (Supplementary Figure 5). These dominant 10 phyla maintained their representation across all lakes sampled, with a variation in relative abundance reflective of the moisture zone. Bacteroidota (0.31 ± 0.08), Proteobacteria (0.24 ± 0.08), and Cyanobacteria (0.19 ± 0.11) were consistently abundant in wet soils. Bacteroidota (0.28 ± 0.09), Proteobacteria (0.18 ± 0.10), and Actinobacteriota (0.15 ± 0.13) were highly abundant in the transition zones. Dry zones were dominated by members of Actinobacteriota (0.30 ± 0.07), Bacteroidota (0.20 ± 0.05), Proteobacteria (0.11 ± 0.03) and Acidobacteriota (0.10 ± 0.06) phyla. In Lake Bonney, taxa affiliated with the phyla Cyanobacteria and Firmicutes presented a relatively high abundance in transition zone (0.16 ± 0.08 and 0.14 ± 0.16, respectively) and dry samples (0.36 ± 0.09 and 0.09 ± 0.09, respectively).

Community beta-diversity from communities from the dry zones of space-for-time transects varied significantly between lakes (PERMANOVA, R² = 0.69, F = 11.55, p-value < 0.01; Supplementary Figure 6A), sharing only 32% of the OTUs (Supplementary Figure 7A). Microbial communities from the wet zones of space-for-time transects were significantly different between lakes (PERMANOVA, R² = 0.73 F = 4.47, p-value < 0.01; Supplementary Figure 6B), sharing only 16% of the OTUs (Supplementary Figure 7B).

Along the space-for-time transects, communities collected in the dry zone were structurally distinct from those collected in the wet zone (PERMANOVA, R² = 0.27, F = 14.11, p-value < 0.01) (Figure 4). For each lake, significant clustering of the communities by zone was identified (Figure 5). However, only the Lake Brownworth Spaulding Pond and Lake Hoare transects met the premise for homogenous dispersion (Figures 5A,B,F; beta-dispersion > 0.05). For all dynamic lake transects (apart from Lake Hoare), the dispersion of communities within the transition zone was significantly different (beta-dispersion < 0.05) (Supplementary Table 2). Except for Lake Hoare, the transition zone community across all dynamic lake transects was segregated into two clusters, with one of the clusters intercepting with the dry soil community cluster on sampling point 3 (ellipses at 95% confidence level) (Figures 5C–E). The latter was not observed in the stable lake transects collected from Lake Brownworth (Figure 5A) or Spaulding Pond (Figure 5B).

Phylogenetic community composition within each sample revealed a trend leaning toward a deterministic assembly of the community governed by environmental selection (NTI > 2) (Supplementary Figure 8A). That effect was higher, and consistent in the dry soils (sampling point 5) and tended to decrease toward the wet zone (sampling point 1) of the gradient, except for Lake Bonney and Lake Vanda, where no trend was identified (Figure 6). In contrast to NTI values, βNTI values distribution showed a median close to 0 but with skewed distribution stretching beyond the + 2 significant threshold (Supplementary Figure 8A). This reflects the influence of deterministic processes for a given pairwise comparison. However, pairwise comparisons regressed against differences in moisture content indicates a large influence of stochasticity with a weak relation to moisture content (Supplementary Figure 8B).

Soil moisture content was the primary non-categorical variable driving microbial community structure patterns (PERMANOVA, R² = 0.25, F = 26.48, p-value < 0.01; Table 2). A random forest approach identified that 64.8% of the variance in the community could be explained by moisture content (RF, no of trees: 10001, mean of squared residuals 19.88). Taxa affiliated to Actinobacteriota, Chloroflexi, Deinococcota, Verrucomicrobiota, Proteobacteria, Acidobacteriota, Abditibacteriota, and Bacteroidota were the major taxa driving differences between wet and dry soils (out-of-bag estimate of error rate for the transect section classification: 24.4%) (Figure 7A). OTUs affiliated to Actinobacteria, order Solirubrobacterales; Chloroflexi, order Kallotenuales; Deinococcota, order Deinococcales; Acidobacteria, order Blastocaellales, and Abditibacteriota, order Abditibacteriales, were exclusively associated to dry and transition zone soils (Figure 7B and Supplementary Table 3). Our model also identified the OTU assigned to Verrucomicrobiota, order Verrucomicrobiales (genus Luteolibacter) as an important OTU to classify wet soils (Figure 7B and Supplementary Table 3). Members of Bacteroidota and Proteobacteria phyla were associated with both dry or wet soil conditions. Bacteroidota taxa related to wet soil conditions belonged to the Flavobacterium and Ferruginobacter genus (order Flavobacteriales and Chitinophagales), whereas the genus Segitobacter (order Chitinophagales) was more abundant in dry soils. Within Proteobacteria, 4 OTUs affiliated to the genus Sphingorhabdus (order Sphingomonadales), Thermomonas and Pseudoxanthomonas (order Xanthomonadales), and genus Brevundimonas (order Caulobacterales) were associated with wet soil conditions, whereas the genus Sphingomonas (order Sphingomonadales) were associated with the dry zone of the gradient (Figure 7B and Supplementary Table 3).

By comparing geochemically defined space-for-time transects from lake systems with relatively stable or dynamic water levels, we were able to resolve the legacy impacts of water presence on microbial community structure (Figure 4). Only within transects defined from Lake Brownworth, and to a lesser extent Spaulding Pond (static transects), did the microbial community display specific structural patterns reflecting the three geochemical zones identified (Figures 5A,B). The high degree of zone specificity (wet, transition, and dry) in community structuring from these lakes aligns with an extended period of more stability due to the presence of an outflow (Levy et al., 2018). In contrast, from the dynamic closed-basin lake transects, the clear spatial segregation of the transition zone communities into one of “dry” or “wet”-like community could reflect a response related to legacy of water level changes in these lakes, not detected by the defined geochemical metrics, which respond on a different time scale. Particularly in Lake Vanda, Lake Bonney, and Lake Fryxell (Figure 5C–E), which have risen between 1.6 and 4.1 m from 2001 to 2017 (Levy et al., 2018), no distinct transitional zone community was recognizable (Figure 7A). Our findings support recent monitoring of Lake Vanda water levels, which indicate that this lake has been continually rising at a rate of 22 cm/year since 1947 (Castendyk et al., 2016). The lack of zone specificity in the microbial community structure at this lake likely reflects that at the time of sampling, microbial communities collected at the geochemically defined transition between transition and dry zones were likely to have been recently wet (potentially hours to days) and remain similar to a dry soil community. In contrast, communities exposed to wetness for more extended periods resemble a more long-term wet soil community. Similar temporal shifts were observed by Ramoneda Massague et al. (2021) in moat samples collected along a depth gradient in Lake Vanda. Microbial assemblages collected from the higher depths diverged from those collected at the lake shore, being likely influenced by pre-existing terrestrial microbes in recently inundated soils (Ramoneda Massague et al., 2021). The consistency of our results across closed-basin lake transects suggests that a geochemically defined space-for-time sampling approach can be applied across polar systems to monitor change through the lasting legacy impacts of water availability on assembly mechanisms that shape extant microbial community structure.

Upon observing different spatial structuring of microbial communities, we evaluated whether microbial community assembly patterns within and between communities along the moisture gradient were driven by deterministic or stochastic processes, or by a combination of both. Contrary to a previous observation (Lee et al., 2018), the NTI values observed for each sample across all 17 geochemically defined space-for-time transects indicate that the degree of phylogenetic relatedness among the individuals in a community is likely driven by environmental selection (NTI > 2) (Figure 6). In the MDV, the dry soils are typically oligotrophic, with minimal water content (<2%) and high salt concentrations (Cary et al., 2010), which collectively exclude a significant range of taxa that cannot tolerate such selective pressures (Chase, 2007; Stegen et al., 2012). Nonetheless, overall patterns of community assembly along increasing differences in moisture content (βNTI values) suggest a considerable but variable influence of stochastic and, to a lesser extent, deterministic processes, with a weak relation to differences in moisture content (Supplementary Figure 8). While this was surprising considering the NTI values, it suggests that, at a regional scale, community assembly is driven by stochastic processes, but phylogenetic differences within communities are driven by deterministic processes that favor the occurrence of closely related taxa. Such a pattern could indicate that microbial traits that confer adaptability are phylogenetically not well-conserved across the region. Alternatively, it could also result from the distribution of spore-forming, or dormant bacteria in dry or highly conductive soils (Lee et al., 2018), or from differences in diversity between lakes or along the gradients. Overall, we suggest that deterministic and stochastic processes play a role in community assembly, however, further extensive analyses are required to understand how deterministic and stochastic components change communities as moisture content changes in the soils.

The impact of moisture as a key driver of microbial diversity was evident at all lakes studied except in Lake Bonney and Lake Vanda (Figure 3A). The presence of water in an available state for life has quantifiable impacts on soil microbial diversity within polar deserts (Niederberger et al., 2015, 2019; Buelow et al., 2016; Bottos et al., 2020). Nonetheless, its positive impact can be suppressed by the effect of soil conductivity (Van Horn et al., 2014; Feeser et al., 2018; Bottos et al., 2020; George et al., 2021). This was the case of Lake Bonney, where the high levels of soil conductivity affected the microbial diversity along the wetness gradient (Figure 3B), creating a niche that only allows for taxa capable of tolerating high osmotic stress (e.g., Firmicutes) to thrive. Lake level rise caused by increased glacial runoff in Lake Bonney facilitates the leaching and movement of accumulated meteoric salts, altering soil geochemistry, particularly pH and conductivity (Barrett et al., 2009). High conductivity measurements in these soils support our evidence for the greater prevalence of deterministic processes driving within community assembly in Lake Bonney and, to a lesser extent, in Lake Vanda, where the moisture content and relatively higher salt concentration in the soils may have a constant selective effect on the community. Both lakes are located further inland, where soils are characteristically older, more developed, and with higher salt concentrations and with different chemistry, compared to coastal soils (Bockheim, 2002). These could also be factors for such disparities in community diversity and composition along transects sampled from these two inland lakes. Although pH has been suggested to be a critical factor driving spatial patterns of microbial communities (Tripathi et al., 2018), in this study, pH changes within and between space-for-time transects did not have a significant effect on community diversity.

Specific phyla and orders of microbial taxa associated with long-term dry and wet conditions can provide a metric to assess early signs of change (Figure 7). Despite the variability in the structure of microbial communities among wet soil samples, OTUs affiliated with Verrucomicrobiota were consistent indicators of long-term wet conditions, whereas OTUs affiliated with Proteobacteria and Bacteroidetes were associated with both wet and dry conditions. On the other hand, Actinobacteria, Acidobacteria, Chloroflexi, and Deinococcota taxa significantly and consistently enriched communities from dry soils, corroborating our previous observations conducted in the MDVs (Niederberger et al., 2015, 2019; Zhang et al., 2020). A recent study by Lee et al. (2018) found similar taxonomic trends, with Acidobacteria and Actinobacteria being less abundant with higher moisture content. We then provide strong evidence that the presence of discrete taxa within the phyla Actinobacteria (family Solirubrobacteraceae), Acidobacteria (subdivision 4 genera Blastocatella), and Deinococcota (genus Truepera) are ubiquitously associated with long-term dry edaphic conditions. This could be a result of metabolic specialization in the shape of very limiting long-term dry, oligotrophic, and ionizing (UV radiation) conditions to life. Taxa belonging to the family Solirubrobacteraceae are currently known for their capability to metabolize atmospheric trace gases contributing to primary production in the dry oligotrophic soils (Ji et al., 2017), and members of the genus Truepera are known for being extremely resistant to ionizing radiation, which is considered a strong stress factor in arid desert soils (Cary et al., 2010). The alleviation of specific stress conditions will impose selective pressures on these highly adapted taxa. As such, any quantifiable taxonomic changes among these taxa can be developed as sentinels for increased temporal exposure to moisture within the MDVs soils, which from a structural standpoint shows a heightened sensitivity of the terrestrial Antarctic microbiome to hydrological changes.

The space-for-time approach employs ecological structures within contemporary environmental gradients to project ecological responses to changes in the environment over time. Primarily applied to study slow ecological processes, this approach is informative when factors that drive community change are equivalent across space and time (Blois et al., 2013). From a microbial ecology perspective, the geochemical profile of a system is often highlighted as an important, if not the primary deterministic force driving community assembly, particularly in ecosystems with low trophic complexity or where diversity is environmentally constrained (Chase, 2007). Yet, sampling designs along spatial gradients are often conducted following distance points (Yang et al., 2014; Niederberger et al., 2015; Yan et al., 2017; Feeser et al., 2018; Lee et al., 2018) rather than being based on the spatial variability of geochemical factors, to which biological communities respond to, along a constrained gradient. The patterns of geochemical drivers across space may not be continuous along a gradient, therefore without characterizing them, geochemical boundaries that may drive community assembly within the gradient are likely to be missed, contributing to extant microbial community heterogeneity.

In an environmentally constrained system, such as the Antarctic polar deserts, we demonstrate that contemporary and historical environmental changes (e.g., water availability) can be resolved by subtle structural and compositional re-arrangements of the microbial communities along geochemically defined environmental gradients. These observations could only be depicted when sampling approaches are determined by the spatial variability of geochemical drivers along environmental gradients and replicated across comparable environments that differ in historical exposure to the driving factor. The validation of our space-for-time transects relies on the observed interdependence between structural and compositional community shifts across space, with the temporal exposure to wet conditions.

Given the accelerated rate of climate change across the globe (IPCC, 2021), the ongoing increasing lake levels (Castendyk et al., 2016; Levy et al., 2018), and the expected increase in snow and ground ice melt (Linhardt et al., 2019) will lead to drastic changes in the MDVs landscape, hydrology and ecology (Fountain et al., 2014). For instance, a recent study conducted in a high Arctic lake demonstrated that the taxonomic and functional diversity of dominant microbes can be impacted by an increase in glacial runoff as a consequence of rapid climate change (Colby et al., 2020). Our study, demonstrates that predictions of a wetter system will directly affect the stability of microbial communities currently adapted to dry and ultra-oligotrophic conditions, resulting in significant compositional and diversity changes across the system. Such environmental and taxonomical changes have the potential to lead to alterations in primary productivity, carbon and nitrogen cycling, and energy fluxes across the system (Coyne et al., 2020; Monteiro et al., 2020; Karen et al., 2021). Whilst microorganisms underpin the stability and functionality of polar deserts ecosystems, these are still rarely considered in policy development and lack any protective status within the Antarctic Treaty (Hughes et al., 2015). This study highlights that microbial communities are highly sensitive to change and closely reflect changes in the environment. Therefore, as in more temperate environments (Cavicchioli et al., 2019), microbial communities should be considered by policy-makers as primary sentinels for current and historical changes in the Antarctic systems.