Unearthing Shifts in Microbial Communities Across a Soil Disturbance Gradient

Permafrost, an important source of soil disturbance, is particularly vulnerable to climate change in Alaska where 85% of the land is underlained with discontinuous permafrost. Boreal forests, home to plants integral to subsistence diets of many Alaska Native communities, are not immune to the effects of climate change. Soil disturbance events, such as permafrost thaw, wildfires, and land use change can influence abiotic conditions, which can then affect active layer soil microbial communities. In a previous study, we found negative effects on boreal plants inoculated with microbes impacted by soil disturbance compared to plants inoculated with microbes from undisturbed soils. Here, we identify key shifts in microbial communities altered by soil disturbance using 16S rRNA gene sequencing and make connections between microbial community changes and previously observed plant growth. Additionally, we identify further community shifts in potential functional mechanisms using long read metagenomics. Across a soil disturbance gradient, microbial communities differ significantly based on the level of soil disturbance. Consistent with the earlier study, the family Acidobacteriaceae, which consists of known plant growth promoters, was abundant in undisturbed soil, but practically absent in most disturbed soil. In contrast, Comamonadaceae, a family with known agricultural pathogens, was overrepresented in most disturbed soil communities compared to undisturbed. Within our metagenomic data, we found that soil disturbance level is associated with differences in microbial community function, including mechanisms potentially involved in plant pathogenicity. These results indicate that a decrease in plant growth can be linked to changes in the microbial community and functional composition driven by soil disturbance and climate change. Together, these results build a genomic understanding of how shifting soil microbiomes may affect plant productivity and ecosystem health as the Arctic warms.

Introduction

Across the sub-arctic and arctic regions, the warming climate is rapidly affecting Alaska’s ecosystems through shifts in disturbance regimes, including increased fire and permafrost thaw (Chapin et al., 2008; Johnstone et al., 2010; Schuur et al., 2015; Chasmer and Hopkinson, 2017). Boreal forests, which represent 90% of the world’s forests, are a complex mosaic of coniferous trees and peatlands (Wolken et al., 2011). Much of the boreal forests across Alaska are underlain by discontinuous permafrost that is not immune to the pressures of climate change. Permafrost thaw is associated with direct and indirect changes in plant communities due to significant shifts in soil hydrology which in turn affect nutrient availability and carbon dynamics (Schuur et al., 2007; Yang et al., 2013; Inglese et al., 2017; Sewell et al., 2020), yet only few studies have explored the biotic mechanism contributing to changes in plant community with thaw (Hewitt et al., 2016; Schütte et al., 2019; Seitz et al., 2021).

Permafrost thaw leads to rapid changes in microbial community composition and function, including shifts in several metabolites and genes involving nitrogen and carbon cycling in response to permafrost thaw, and during a thaw event, community function within permafrost quickly converges to that of the active layer (Mackelprang et al., 2011; Coolen and Orsi, 2015; Monteux et al., 2018; Johnston et al., 2019; Messan et al., 2020; Saidi-Mehrabad et al., 2020). Consistent with the pattern of permafrost thaw affecting soil microbial communities, Inglese et al. (2017) found that active layer detachments, a form of permafrost disturbance, significantly affect fungal and Archaeal community composition of the active layer. They observed greater proportions of Nitrososphaerales, an ammonia oxidizing Archaea, in disturbed communities compared to undisturbed, again indicating a shift in nitrogen cycling following permafrost thaw. Within fungal communities, Inglese et al. (2017) saw a large decrease in ericoid mycorrhizal fungi and Ascomycota overall in disturbed soils compared to undisturbed communities. These groups are often found associated with arctic ericoid plant species. These authors suggest that both the increase of ammonia-oxidizing Archaea and the reduction in hyphal fungi could exacerbate further permafrost thaw and landscape change in the arctic.

In a previous study, we found that increased soil disturbance was associated with decreased plant productivity across multiple key boreal plant species, and that this decrease may be mediated by underlying differences in soil microbial community composition (Seitz et al., 2021). The microbial community analysis in Seitz et al. (2021) was based on long read sequence analysis using Oxford Nanopore sequencing technology. However, the usefulness of long read sequence data is not yet well established for environmental studies. Here, we had the unique opportunity to determine whether 16S rRNA gene sequencing would show consistent differences in community composition compared to the previously obtained nanopore sequences, as well as increase our understanding of how disturbance affects alpha and beta diversity metrics. In this study, we combined microbial community composition and functional analysis, and linked belowground changes with changes in plant productivity, which statistically results in a comprehensive understanding of disturbance-induced changes in belowground processes affecting plant communities.

In this study, we determined differences in microbial community composition and function across the disturbance gradient using both 16S rRNA gene sequencing and long read metagenomic sequencing based on nanopore sequencing technology. Further, we used linear regression analysis to determine whether within disturbance level variation can predict plant productivity previously observed in the plant growth experiment (Seitz et al., 2021). We hypothesized that both alpha and beta diversity would increase with disturbance, and furthermore, not only would microbial community composition differ across disturbance levels, but also potential community function. The data of this study significantly add to our previous study (Seitz et al., 2021) providing a more comprehensive understanding of how disturbances affect belowground processes and that these changes in microbial community composition and function indeed significantly impact plant communities.

Materials and Methods

Sample Site Description

The soil samples we collected and used for this study were previously described in Seitz et al., 2021. We collected samples from interior Alaska from the Fairbanks Permafrost Experiment Station (FPES; 64.877°N, 147.670°W). This forested site was established and disturbed in 1946 as part of the United States Army Corps of Engineers Cold Regions Research and Engineering Laboratory, which aimed to simulate and study how disturbance events influence the ecology of the boreal forest. FPES contains three, 3,721 m² Linell plots, described in Linell, 1973 of varying degrees of soil disturbance. When FPES was established, the first plot, referenced as undisturbed (UD), was left untouched. The second plot (semi-disturbed, SD) was cleared of all above ground vegetation and trees, while roots and organic soil layers stayed undamaged. The third plot (most disturbed, MD) was stripped entirely of surface level vegetation and organic material until the mineral soil was reached. Since the original disturbance event, permafrost in the SD plot degraded to 4.7 m below the surface soils, and the permafrost in the MD plot degraded down to 9.8 m below the surface (Douglas et al., 2008).

Fairbanks Permafrost Experiment Station is an example of a Taiga boreal forest found in subarctic regions across the world. The untouched UD plot is characterized by a thick covering of black spruce (Picea Mariana) with white spruce (Picea glauca) scattered throughout. The UD understory is largely dominated by low-bush cranberry (lingonberry, Vaccinium vitis-idaea) and Labrador tea (Rhododendron groendlandicum), and the ground is continuously covered with Sphagnum and feather moss, and lichen. The SD plot is characterized by a highly mixed stand of black spruce, white spruce, Alaskan birch (Betula neoalaskana), and willow (Salix spp.). SD trees are taller than those in the MD plot which is characterized by young willow, black spruce, and Alaskan birch (Douglas et al., 2008).

Sample Collection and DNA Extraction

To assess how soil microbial communities are shifting based on previous soil disturbance, we collected a total of 48 soil samples on 28 May 2018, with 16 each from UD, SD, and MD plots. The collection of these cores and their original use was described in Seitz et al. (2021). In brief, we extracted cores from a selection of quadrats throughout each disturbance plot to sample the within disturbance level heterogeneity. Using sterile technique, we removed the top layer of moss at each sampling point and then collected the top layer of soil using a 4.5 cm diameter by 10 cm height soil corer. We stored samples in a cooler throughout the duration of sample collection before they were transported back to the lab. After the samples arrived at the lab, we mechanically homogenized each individual soil core and stored them at 4°C.

After homogenizing each soil sample, we extracted DNA from 0.25 g of soil using the Qiagen PowerSoil Kit (Qiagen Inc., United States) according to the manufacturer’s instructions. We then checked the quality and concentration of DNA using a NanoDrop One spectrophotometer (ThermoFisher Science, United States) and a Qubit 4.0 fluorometer (Fisher Science, United States).

PCR and 16S rRNA Gene Sequencing and Analyses

We determined microbial community composition based on 16S rRNA gene sequencing, which allowed us to determine taxonomy at a finer resolution compared to our previous study and gave us the unique opportunity to compare the commonly used standard of 16S rRNA gene sequences to the not as well-established long read based on nanopore sequencing (Seitz et al., 2021). After extracting DNA, we performed a 1:10 dilution of all samples to be used for PCR. We amplified the V4 region using dual-indexed 515F (Parada et al., 2016) and 806R (Apprill et al., 2015) primers following the EMP PCR protocol¹ to be used for 16S rRNA gene amplicon analysis. Following the first round of amplification, we used gel electrophoresis to determine the presence of PCR products. We ran duplicate PCR reactions and then pooled the two reactions for each unique sample. The amplicons were then sequenced on an Illumina MiSeq using v3 reagents at the Institute of Arctic Biology Genomics Core Lab.

Upon obtaining the 16S rRNA amplicon data, we first demultiplexed all reads using Mr_Demuxy (version 1.2.0, https://pypi.org/project/Mr_Demuxy/1.2.0/). We then processed and analyzed the paired-end, demultiplexed reads using the Quantitative Insights into Microbial Ecology (QIIME2) framework (version 2020.8.0). We used the QIIME 2 plugin DADA2² to obtain a table of representative sequences and frequencies, remove low quality regions of reads, and to merge our paired reads for further processing. Next, we constructed a rooted phylogenetic tree from representative sequences using the QIIME 2 alignment plug in and then calculated microbial diversity metrics based on amplicon sequence variants (ASVs) allowing taxonomic classification at a finer resolution compared to the long read data based on nanopore sequencing used in our previous study (Seitz et al.). After constructing the tree, we assigned taxonomy with the QIIME 2 taxonomy classifier plugin (Quast et al., 2012; Bokulich et al., 2018) using a Naïve Bayes classifier pre-trained on the Silva 138 99% OTUS (515F/806R region). We next filtered out any chimeric, mitochondrial, and chloroplast reads from our samples using the QIIME 2 plugin feature table.

Amplicon Statistical Analyses

We determined alpha and beta diversity, tested whether diversity differed significantly across disturbance levels, and determined whether microbial community composition was significantly associated with changes in plant productivity. This type of analysis was absent in our previous study (Seitz et al., 2021), and 16S rRNA gene sequences enabled us to complete these analyses based on finer taxonomic resolution (ASV). All statistical analyses performed in R utilized version 3.6.1. After obtaining our filtered feature table, we used qiime2R (version 0.99.35; Bisanz, 2018) to import our QIIME 2 data and artifacts into R Studio. We converted our taxonomic data into a phyloseq (version 1.30.0; McMurdie and Holmes, 2013) object and filtered out one sample with low sampling depth (less than 14,377 reads). We transformed our data into relative abundance measures. We used vegan (version 2.5.7; Oksanen et al., 2020) to calculate alpha diversity metrics (Pielou’s evenness and Faith’s phylogenetic diversity) based on ASVs. We visualized our results with ggplot2 (version 3.3.3; Wickham, 2016).

To test the significance of overall differences between UD, SD, and MD soil communities, we calculated the Bray–Curtis similarity index to compare community beta diversity using amplicon sequence variants (ASVs). We tested the effects of disturbance on community composition using the Permutational Multivariate ANOVA (PERMANOVA) using the “adonis” function in vegan (10,000 permutations). We visualized beta diversity differences using nonmetric multidimensional scaling (NMDS) and ggplot2.

To visualize community membership across soil samples and taxonomic levels, we generated heatmaps using the R package gplots (version 3.1.1; Warnes et al., 2020). At the phylum level, we filtered out any phyla that were not present in samples at higher than 1% relative abundance. When visualizing at the family level, we filtered out any family that was present at a maximum abundance of less than 5% across all samples.

To examine the statistical relationship between microbial community composition and resulting plant productivity Seitz et al. (2021) observed in their growth experiment, we utilized linear regression with the lm function in R (R Core Team, 2019). Since the first NMDS axis from our 16S rRNA amplicon analysis describes microbial community variation across soil disturbance level, we extracted the scores from the axis of variation, NMDS1, to test whether community variation within disturbance level explains leaf count and plant height from Seitz et al. (2021). For example, for bog blueberry plants, we performed linear regression using MD community variation NMDS scores and plant height as the growth measure. We repeated the regression using SD and then UD scores. We conducted this series of analyses for each plant type and growth measure (height and leaf count). We then adjusted p values using the Benjamini-Hochberg method with a false discovery rate of 0.25 on the basis that follow up studies are relatively low cost (Benjamini and Hochberg, 1995). We then used ggplot2 (version 3.3.3) to visualize significant relationships between plant productivity and microbial community variation within each soil disturbance level.

Next, we compared community composition estimates between the gold standard 16S rRNA amplicon sequencing and long read nanopore sequencing using the base R package, stats (version 3.6.1; R Core Team, 2019). We calculated the Pearson’s correlation coefficient estimates between 16S and long read relative abundances for 10 specific taxonomic groups previously highlighted as biomarkers of either UD or MD communities in Seitz et al. (2021).

Functional Analysis

We performed community functional analyses on previously collected metagenomic data from our soil cores described in Seitz et al., 2021 (ENA Project: PRJEB42020). Determining microbial function based on the metagenomic data was an aspect lacking in our previous study (Seitz et al., 2021), but essential to obtain a more complete understanding on the effects of disturbances such as permafrost thaw on microbial communities. We utilized the MEGAN-LR (long read) pipeline (Huson et al., 2018) following the parameters described by Arumugam et al. (2019). MEGAN-LR is a tool that performs both taxonomic and functional community classification. Briefly, this pipeline aligns long reads against an NCBI-nr database using DIAMOND (Buchfink et al., 2015) by performing a frame-shift aware DNA-to-protein alignment. The alignments are then processed by MEGAN using an LCA-based algorithm for taxonomic binning and an interval-tree based algorithm for functional binning.

Following classification of potential functions via MEGAN-LR, we utilized the MEGAN comparison function to normalize read counts and compare functions across samples. We then extracted the evolutionary genealogy of genes: Non-supervised Orthologous Groups (EGGNOG)-Clusters of Orthologous Group (COG) functions into a matrix. We calculated differences in functional diversity metrics using the vegan R package (version 2.5.7; Oksanen et al., 2020) and visualized differences in function based on the identified COG functions using nonmetric dimensional scaling (NMDS) using ggplot (version 3.3.3; Wickham, 2016). We tested the effects of disturbance level through a PERMANOVA using the “adonis” function of vegan. To identify potential functional indicators of each disturbance level (UD, SD, and MD), we used the R package indicspecies (version 1.7.9; De Caceres and Legendre, 2009). In this application, indicator species analysis is using indices of a gene’s or taxon’s (in this case gene) relative abundance and their occurrence to estimate the strength of its association with specific groups (in this instance soil disturbance level). For the indicator analysis, we used only the identified COGs that were identified via MEGAN-LR. Within indicspecies, we used the function “multipatt” to determine a list of COGs that are significantly associated with each disturbance level. We specified 999 random permutations. We visualized the results using the R package gplots (version 3.1.1; Warnes et al., 2020).

Results

Microbial Community Composition Differs Based on Soil Disturbance Level

We identified a total of 17,005 ASVs within the 47 retained samples with a feature count ranging from 14,076 to 66,994, and a mean feature count of 50,943. After filtering out ASVs that were not seen more than three times in at least 20% of samples, we detected 1,005 ASVs. Using the filtered ASV data, we compared both alpha and beta diversity to test for differences in bacterial community across soil disturbance level. We found alpha diversity increased significantly within communities as the level of soil disturbance increased (Figure 1). A Kruskal-Wallis test revealed that the Faith’s Phylogenetic Diversity within MD soil communities was significantly higher than in either SD or UD communities (Table 1). We also found that Pielou’s Evenness did not differ significantly between MD and SD communities but did increase from UD to MD and SD soils.

FIGURE 1

Figure 1. Boxplots of (A) Faith’s Phylogenetic Diversity and (B) Pielou’s Evenness.

TABLE 1

Table 1. All groups and pairwise comparison of alpha diversity by Kruskal-Wallis.

We found that microbial community composition differed significantly at the ASV level across the three soil disturbance levels (Figure 2; stress = 0.055; PERMANOVA F_2,44 = 28.682, p < 0.001). Communities clustered based on disturbance level along the NMDS first axis. MD communities clustered distinctly from both SD and UD communities, while SD and UD communities displayed some overlap. The variation in microbial community composition along axis 2 was not associated with separation of the microbial communities across the disturbance level groups. To compare the taxonomic results of 16S and metagenomic sequencing, we completed correlation analyses (Pearson’s correlation coefficient) on the relative abundances of biomarkers highlighted in Seitz et al. (2021). We found highly significant positive correlations between 16S and metagenomic relative abundances of eight biomarkers ranging from 0.641 to 0.942 with an average of 0.808 (Table 2). Two biomarkers, the phylum Proteobacteria and the order Myxococcales, did not show significant correlations between sequencing methods.

FIGURE 2

Figure 2. Nonmetric multidimensional scaling (NMDS) based on Bray–Curtis dissimilarity distances showing the differences in beta diversity between undisturbed (UD), semi disturbed (SD), and most disturbed (MD) soil communities. Points colored by the level of Fairbanks Permafrost Experiment Station (FPES) soil disturbance with blue = UD (n = 16), gold = SD (n = 15), and red = MD (n = 16).

TABLE 2

Table 2. Pearson’s correlation estimates (r) between 16S and long read relative abundances.

Across all three disturbance levels, we observed shifts in community membership and abundances. We identified 11 phyla that were present at maximum abundances greater than 5% within all samples. Out of those 11 phyla, Proteobacteria were the most abundant phylum, present at mean relative abundances of 40.7% (range 25.2%–58.8%) in UD, 36.1% (range 27.5%–43.8%) in SD, and 38.0% (range 29.9%–46.4%) in MD soil communities. Looking further at the taxonomic composition of the soil communities, we found 17 groups classified to the family level that were present at maximum abundances greater than 5% (Figure 3). Of those, Acidobacteriaceae was the most abundant family in UD, present at a mean relative abundance of 16.5% (range 5.46%–24.7%) compared to a mean of 0.33% (range 0.00%–3.21%) in MD communities. The most abundant families found in MD communities were Comamonadaceae and Nitrosomonadaceae, both Proteobacteria, present at mean relative abundances of 3.33% (range 2.0%–5.7%) and 6.53% (range 3.8%–10.3%), respectively. Within UD communities, Comamonadaceae and Nitrosomonadaceae were present at mean relative abundances of 0.28% (range 0.03%–1.57%) and 0.05% (range 0.00%–0.51%), respectively. Within the phylum Acidobacteria, Subgroup 2 and “uncultured” bacteria, and WD260, a group of uncultured bacteria belonging to the phylum Proteobacteria were present at relative abundances greater than 5% throughout our samples.

FIGURE 3

Figure 3. Heatmap of bacterial families present at maximum abundances greater than 5%. Each row corresponds to a family, and each column corresponds to an individual soil core community. The top row signifies the level of soil disturbance of each core with blue = undisturbed (UD; n = 16), gold = semi disturbed (SD; n = 15), and red = most disturbed (MD; n = 16). The rows are clustered by Bray-Curtis distance and the color of each box displays the relative abundance of the taxon within that core.

Microbial Community Composition as a Predictor for Plant Growth

Following observing significant changes in microbial community diversity depending on the level of soil disturbance, we sought to determine whether variation within disturbance level community composition could predict the growth of plants from an earlier experiment described in Seitz et al. (2021). Within each disturbance level we regressed the productivity measures (height or leaf count) against the NMDS1 values for each corresponding microbial community used as an inoculant in the plant growth experiment (Seitz et al., 2021). We observed significant negative linear relationships between MD microbial community composition and plant growth within bog blueberry height, low-bush cranberry height and leaf count, and Labrador tea height and leaf count (Table 3; Figure 4). For bog blueberry, low-bush cranberry, and Labrador tea plants grown with MD microbes, plant growth continued to decrease as microbial communities became more extreme and more positive along the NMDS1 axis (NMDS1 > 0.5. Microbial community composition within SD and UD communities was not a significant predictor of height or leaf count in any plant species except for black spruce (height) from Seitz et al., 2021; Table 3).

TABLE 3

Table 3. Results from linear models of plant growth measures predicted by 16S NMDS1.

FIGURE 4

Figure 4. Scatterplots showing the effects of microbial community variation and plant productivity (16S rRNA amplicon NMDS1 and height/leaf count) of bog blueberry (A,B), low-bush cranberry (C,D), and Labrador tea (E,F). Each point represents a soil community, colored by the level of FPES soil disturbance with blue = undisturbed (UD; n = 16), gold = semi disturbed (SD; n = 15), and red = most disturbed (MD; n = 16). Solid lines represent a significant linear relationship.

Differences in Community Function Based on Soil Disturbance Level

Based on the observed EggNOG annotations from our metagenomic data, we identified 2,891 clusters of orthologous genes (COGs) and 7,980 non-supervised orthologous groups across all samples. Using Bray-Curtis distances and a PERMANOVA analysis, we found that potential community function varied significantly based on the level of soil disturbance (Figure 5; stress = 0.14; PERMANOVA F_2,44 = 0.17431, p = 0.001).

FIGURE 5

Figure 5. NMDS of Bray–Curtis dissimilarity distances of community function. Functions based on MEGAN-LR evolutionary genealogy of genes: Non-supervised Orthologous Groups (EGGNOG) Clusters of Orthologous Group (COG) classifications. Points colored by the level of FPES soil disturbance with blue = undisturbed (UD; n = 16), gold = semi disturbed (SD; n = 15), and red = most disturbed (MD; n = 16).

To identify specific functions that are characteristic of the three levels of soil disturbance, we performed an indicator analysis on the observed potential genomic functions of our microbial communities. After filtering the results for a high specificity (A ≥ 0.8) and sensitivity (B ≥ 0.5), we found 14 COG indicators, with 11 belonging to MD, two SD, and one UD (Supplementary Figure S1; Table 4). Within the indicator COGs, we identified a variety of clusters. The only COG identified from UD encodes for a chloramphenicol acetyltransferase an enzyme involved in the detoxification of the antibiotic chloramphenicol associated with resistance to chloramphenicol in bacteria. The COGs from SD were associated with basal transcription factors and nucleic acid binding. From MD, we identified a variety of COGs including those encoding for a for alpha-N-arabinofurandase (ABF), a cell-wall degrading enzyme, disulfide reductases, and a hydrogenase iron–sulfur subunit involved in methanogenesis and a nitrate reductase beta subunit necessary in denitrification.

TABLE 4

Table 4. Highly sensitive and specific indicator functions from EggNOG annotations.

Discussion

This study provides a thorough microbial analysis including determining differences in microbial community composition at a fine taxonomic resolution, in alpha diversity, and in microbial function based on indicator species analysis results. Combined with the findings in our previous study (Seitz et al., 2021), our data give a more comprehensive picture of how disturbance affects microbial communities and mediates changes in plant communities.

Supporting our hypothesis that soil disturbance drives community shifts; we found that microbial community composition significantly differed with level of soil disturbance at the ASV level (Figures 1, 2) despite within site heterogeneity especially in UD and SD. This is consistent with our previous findings that community composition differed across these three disturbance levels (Seitz et al., 2021). We continued to look at these shifts in taxonomic composition through the lens of biomarkers previously identified in our metagenomic data as either over- or underrepresented in UD and MD communities (Seitz et al., 2021). We found that the relative abundances of these biomarkers were significantly positively correlated between our 16S and metagenomics datasets, further underlining that long read metagenomic sequences are a useful tool to assess differences in microbial community composition across different environmental conditions (Seitz et al., 2021).

We detected significant increases in alpha diversity with increasing in soil disturbance. This is a phenomenon that is observed within many ecosystems subject to many different types of disturbance, suggesting that disturbance promotes diverse communities due to the creation of new niches within the environment (Violle et al., 2010; Galand et al., 2016). Our findings of microbial community diversity increasing with disturbance are consistent with long-term (Feng et al., 2020) and short-term (Luláková et al., 2019) studies of the effects of warming soils and permafrost thaw on active layer bacterial communities. Feng et al. (2020) and Luláková et al. (2019) both saw alpha diversity measures significantly increase as topsoil temperatures increased. This trend of alpha diversity increasing along with soil disturbance has also been demonstrated in soils from a variety of land and forest types that have been subject to physical disturbances, such as logging, agricultural conversion, and fires (Chatterjee et al., 2008; de Carvalho et al., 2016; Shen et al., 2016; Zhou et al., 2018). The physical soil disturbance and resulting deep permafrost thaw that the MD level plot underwent has created a heterogenous landscape (Douglas et al., 2008), which may be contributing to this increased diversity by affecting the ecological niches and competition ability of bacteria for resources, similar to what can occur after forest fires (Shen et al., 2016).

We observed a lower community evenness that within our UD communities compared to our MD and SD communities. The relative abundance of Acidobacteria decreased drastically while relative abundance of Nitrosomanadaceae and Comamonadaceae increased with increasing soil disturbance. When ASVs were grouped at the family level, Acidobacteriaceae was the most abundant bacterial family within the UD communities, compared to being nearly absent in MD communities. Members of Acidobacteriaceae are typically found to be present at higher abundances in acidic soils that are carbon poor (Kielak et al., 2016; Osburn et al., 2019). In contrast, Nitrosomonadaceae, a copiotrophic bacterial family of ammonia oxidizers, was found to be the most abundant family within the MD soil communities, followed by Comamonadaceae. Nitrosomonadaceae are ammonia oxidizing bacteria (AOB) that are commonly found in carbon-rich environments such as wastewater and soils (Prosser et al., 2014). An increase of copiotrophic bacteria with disturbance is consistent with the body of literature (Masse et al., 2017; Mickan et al., 2019; Osburn et al., 2019). Masse et al. (2017) previously detected that oligotrophic bacteria, such as Acidobacteriaceae, were more abundant in natural boreal forest soils, compared to soils that were reclaimed or reconstructed after disturbance which displayed higher abundances of copiotrophic bacteria.

We found that the microbial community composition across soil disturbance is a significant predictor of plant productivity measures for bog blueberry (Vaccinium uliginosum), low-bush cranberry (Vaccinium vitis-ideae), and Labrador Tea (Rhododenron groenlandicum). Specifically, we can predict that plants inoculated with a microbial community indicative of the most extreme disturbance, will demonstrate decreased growth compared to plants inoculated with communities from less disturbed soils. This significant relationship between microbial community composition and plant productivity was observed within plants inoculated with MD microbes, whereas we saw no significant relationships between plant productivity of plants inoculated with microbial communities from SD or UD and microbial community composition. This suggests that there is a threshold effect occurring across the disturbance levels, and that soil microbial communities subject to disturbance resulting in permafrost thaw may not have the capacity to influence plant productivity of these boreal species until affected by an event of large enough magnitude. This microbial threshold we are observing could potentially be linked to a disturbance threshold of deep permafrost thaw. Over the past 20 years, researchers have observed accelerated permafrost thaw within the northern boreal Taiga plains likely due to a combination of hydro-climactic changes and increases in annual air temperatures (Jorgenson et al., 2006; Lara et al., 2016; Chasmer and Hopkinson, 2017). As permafrost thaw is accelerated and other disturbance events are increasing in the artic, we may start to observe the effects through changes in soil microbial communities and plant productivity.

Our previous study (Seitz et al., 2021) lacked any analysis of community function, and based on the disturbance-associated community membership shifts we observed, we predicted that observed functional community will shift as well. Using the metagenomic sequencing data, we identified functional indicators across all our FPES sites representative of each level of soil disturbance, with most indicators belonging to MD soils. The gene cluster with the highest specificity within MD is COG3940, which encodes for alpha-N-arabinofurandase (ABF). ABFs are mainly extracellular enzymes that degrade lignocelluloses and cell walls (De Ioannes, 2000; Numan and Bhosle, 2006), and have been shown to play a role in triggering plant immune responses. In 2016, Wu et al. found that a novel ABF plays a critical role in the pathogenicity of the fungal pathogen, Magnaporthe oryzae, that causes rice blast disease. The secreted ABF degrades cell walls, further leading to a decrease in productivity in infected plants (Wu et al., 2016). The presence of a gene cluster encoding for alpha-N-arabinofurandase in our MD communities provides a possible mechanism behind the decrease in plant productivity demonstrated in Seitz et al. (2021), yet further studies involving gene expression and enzyme activity would be required to identify if this is part of a mechanism at play.

The critical gene involved in ammonia oxidization (amoA; Prosser et al., 2014; Bárta et al., 2017) was present across all our sample communities, however, this gene cluster was not identified as a functional indicator of MD communities where Nitrosomonadaceae were very abundant. Curiously, we identified COG1140, a gene cluster encoding for a nitrate reductase beta subunit indicative of denitrification, to be a highly specific and significant indicator of MD soil communities. Nitrate losses within soils have long been known to follow disturbance events in forests such as fire, disease, or clear cutting (Vitousek et al., 1979; Neary et al., 1999; Hobara et al., 2001; Pardo et al., 2002). However, comparatively little is known about denitrification in the boreal system. This potential loss of available nitrate could be part of the mechanism at play in the reduced plant productivity that Seitz et al. (2021) observed within plants grown in MD community inoculated soils.

As climate change and anthropogenic-driven soil disturbance events including wildfires (Calef et al., 2008; Chapin et al., 2008; Partain et al., 2016) and permafrost thaw (Jorgenson et al., 2006; Lara et al., 2016; Chasmer and Hopkinson, 2017) are increasing across the north, we need to continue to study how disturbance influences the microbiomes of active layer soils that are both directly and indirectly involved in mediating plant structure within boreal forests. In this study, we show that active layer soil microbial communities are influenced by the initial disturbance event and permafrost thaw at FPES, however we cannot conclude whether specific environmental factors are associated with community shifts. We observed significant shifts in microbial community composition, in microbial diversity, as well as microbial community function. Our evidence for microbial community composition being a significant predictor of decreased plant productivity in response to the most severe soil disturbance associated with significant permafrost thaw suggests that this effect may occur once a threshold level of disturbance is reached. These consistent effects of disturbance on active layer microbial communities, combined with the knowledge of how they in turn can significantly affect plant productivity (Seitz et al., 2021), highlight the need for more research into how on-going and future disturbance events will influence the ecological dynamics in boreal forests.

Data Availability Statement

The original contributions presented in the study are publicly available. This data can be found at: PRJEB42020.

Author Contributions

US and DD contributed to the conception, design, and analysis of the study. TS performed the data collection, statistical analysis and wrote the first draft of the manuscript. All authors contributed to the article and approved the submitted version.

Funding

Research reported in this publication was supported by an Institutional Development Award (IDeA) from the National Institute of General Medical Sciences of the National Institutes of Health under grant number P20GM103395 as well as under three linked awards numbers RL5GM118990, TL4GM118992, and 1UL1GM118991.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We would like to thank Anne-Lise Ducluzeau and Tracie Haan for laboratory support, as well as Scout McDougall and Jennie Humphrey for support with data collection. We would also like to thank Christa Mulder and Brandon Briggs for their feedback and guidance on analysis and writing. The authors thank to Tom Douglas from the Cold Regions Research and Engineering Laboratory (CRREL)-Alaska for providing access to the field site, and many thanks to the reviewers who provided helpful feedback on earlier versions of this work. This research was supported by the Alaska BLaST program, the Institute of Arctic Biology, and Alaska INBRE.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2022.781051/full#supplementary-material

Footnotes

1. ^https://earthmicrobiome.org/protocols-and-standards/16s/

2. ^https://benjjneb.github.io/dada2/index.html

References

Apprill, A., McNally, S., Parsons, R., and Weber, L. (2015). Minor revision to V4 region SSU rRNA 806R gene primer greatly increases detection of SAR11 bacterioplankton. Aquat. Microb. Ecol. 75, 129–137. doi: 10.3354/ame01753

CrossRef Full Text | Google Scholar

Arumugam, K., Bağcı, C., Bessarab, I., Beier, S., Buchfink, B., Górska, A., et al. (2019). Annotated bacterial chromosomes from frame-shift-corrected long-read metagenomic data. Microbiome 7:61. doi: 10.1186/s40168-019-0665-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Bárta, J., Tahovská, K., Šantrůčková, H., and Oulehle, F. (2017). Microbial communities with distinct denitrification potential in spruce and beech soils differing in nitrate leaching. Sci. Rep. 7:9738. doi: 10.1038/s41598-017-08554-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Benjamini, Y., and Hochberg, Y. (1995). Controlling the false discovery rate: a practical and powerful approach to multiple testing. J. R. Stat. Soc. Ser. B Methodol. 57, 289–300. doi: 10.1111/j.2517-6161.1995.tb02031.x

CrossRef Full Text | Google Scholar

Bisanz, J. E. (2018). Qiime2R: Importing QIIME2 Artifacts and Associated Data into R Sessions. Available at: https://github.com/jbisanz/qiime2R (Accessed March 16, 2022).

Google Scholar

Bokulich, N. A., Kaehler, B. D., Rideout, J. R., Dillon, M., Bolyen, E., Knight, R., et al. (2018). Optimizing taxonomic classification of marker-gene amplicon sequences with QIIME 2’s q2-feature-classifier plugin. Microbiome 6:90. doi: 10.1186/s40168-018-0470-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Buchfink, B., Xie, C., and Huson, D. H. (2015). Fast and sensitive protein alignment using DIAMOND. Nat. Methods 12, 59–60. doi: 10.1038/nmeth.3176

PubMed Abstract | CrossRef Full Text | Google Scholar

Calef, M. P., McGuire, A. D., and Chapin, F. S. (2008). Human influences on wildfire in Alaska from 1988 through 2005: an analysis of the spatial patterns of human impacts. Earth Interact. 12, 1–17. doi: 10.1175/2007EI220.1

CrossRef Full Text | Google Scholar

Chapin, F. S., Trainor, S. F., Huntington, O., Lovecraft, A. L., Zavaleta, E., Natcher, D. C., et al. (2008). Increasing wildfire in Alaska’s boreal Forest: pathways to potential solutions of a wicked problem. Bioscience 58, 531–540. doi: 10.1641/B580609

CrossRef Full Text | Google Scholar

Chasmer, L., and Hopkinson, C. (2017). Threshold loss of discontinuous permafrost and landscape evolution. Glob. Chang. Biol. 23, 2672–2686. doi: 10.1111/gcb.13537

PubMed Abstract | CrossRef Full Text | Google Scholar

Chatterjee, A., Vance, G. F., Pendall, E., and Stahl, P. D. (2008). Timber harvesting alters soil carbon mineralization and microbial community structure in coniferous forests. Soil Biol. Biochem. 40, 1901–1907. doi: 10.1016/j.soilbio.2008.03.018

CrossRef Full Text | Google Scholar

Coolen, M. J. L., and Orsi, W. D. (2015). The transcriptional response of microbial communities in thawing Alaskan permafrost soils. Front. Microbiol. 6:197. doi: 10.3389/fmicb.2015.00197

PubMed Abstract | CrossRef Full Text | Google Scholar

De Caceres, M., and Legendre, P. (2009). Associations between species and groups of sites: indices and statistical inference. Ecology 90, 3566–3574. doi: 10.1890/08-1823.1

PubMed Abstract | CrossRef Full Text | Google Scholar

de Carvalho, T. S., Jesus, E. D., Barlow, J., Gardner, T. A., Soares, I. C., Tiedje, J. M., et al. (2016). Land use intensification in the humid tropics increased both alpha and beta diversity of soil bacteria. Ecology 97, 2760–2771. doi: 10.1002/ecy.1513

PubMed Abstract | CrossRef Full Text | Google Scholar

De Ioannes, P. (2000). An α-L-arabinofuranosidase from Penicillium purpurogenum: production, purification and properties. J. Biotechnol. 76, 253–258. doi: 10.1016/S0168-1656(99)00190-X

PubMed Abstract | CrossRef Full Text | Google Scholar

Douglas, T., Kanevskiy, M., Romanovsky, V., Shur, Y., and Yoshikawa, K. (2008). “Permafrost Dynamics at the Fairbanks Permafrost Experimental Station Near Fairbanks, Alaska.” in Ninth International Conference on Permafrost. 29 June–3 July, 2008; 373–378.

Google Scholar

Feng, J., Wang, C., Lei, J., Yang, Y., Yan, Q., Zhou, X., et al. (2020). Warming-induced permafrost thaw exacerbates tundra soil carbon decomposition mediated by microbial community. Microbiome 8:3. doi: 10.1186/s40168-019-0778-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Galand, P. E., Lucas, S., Fagervold, S. K., Peru, E., Pruski, A. M., Vétion, G., et al. (2016). Disturbance increases microbial community diversity and production in marine sediments. Front. Microbiol. 7:1950. doi: 10.3389/fmicb.2016.01950

PubMed Abstract | CrossRef Full Text | Google Scholar

Hewitt, R. E., Bennett, A. P., Breen, A. L., Hollingsworth, T. N., Taylor, D. L., Chapin, F. S., et al. (2016). Getting to the root of the matter: landscape implications of plant-fungal interactions for tree migration in Alaska. Landsc. Ecol. 31, 895–911. doi: 10.1007/s10980-015-0306-1

CrossRef Full Text | Google Scholar

Hobara, S., Tokuchi, N., Ohte, N., Koba, K., Katsuyama, M., Kim, S. J., et al. (2001). Mechanism of nitrate loss from a forested catchment following a small-scale, natural disturbance. Can. J. For. Res. 31, 1326–1335. doi: 10.1139/x01-064

CrossRef Full Text | Google Scholar

Huson, D. H., Albrecht, B., Bağcı, C., Bessarab, I., Górska, A., Jolic, D., et al. (2018). MEGAN-LR: new algorithms allow accurate binning and easy interactive exploration of metagenomic long reads and contigs. Biol. Direct 13:6. doi: 10.1186/s13062-018-0208-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Inglese, C. N., Christiansen, C. T., Lamhonwah, D., Moniz, K., Montross, S. N., Lamoureux, S., et al. (2017). Examination of soil microbial communities after permafrost thaw subsequent to an active layer detachment in the high Arctic. Arct. Antarct. Alp. Res. 49, 455–472. doi: 10.1657/AAAR0016-066

CrossRef Full Text | Google Scholar

Johnston, E. R., Hatt, J. K., He, Z., Wu, L., Guo, X., Luo, Y., et al. (2019). Responses of tundra soil microbial communities to half a decade of experimental warming at two critical depths. PNAS 116, 15096–15105. doi: 10.1073/pnas.1901307116

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnstone, J. F., Chapin, F. S., Hollingsworth, T. N., Mack, M. C., Romanovsky, V., and Turetsky, M. (2010). Fire, climate change, and forest resilience in interior Alaska, this article is one of a selection of papers from The dynamics of change in Alaska’s boreal forests: resilience and vulnerability in response to climate warming. Can. J. For. Res. 40, 1302–1312. doi: 10.1139/X10-061

CrossRef Full Text | Google Scholar

Jorgenson, M. T., Shur, Y. L., and Pullman, E. R. (2006). Abrupt increase in permafrost degradation in Arctic Alaska. Geophys. Res. Lett. 33:L02503. doi: 10.1029/2005GL024960

CrossRef Full Text | Google Scholar

Kielak, A. M., Barreto, C. C., Kowalchuk, G. A., van Veen, J. A., and Kuramae, E. E. (2016). The ecology of Acidobacteria: moving beyond genes and genomes. Front. Microbiol. 7:744. doi: 10.3389/fmicb.2016.00744

PubMed Abstract | CrossRef Full Text | Google Scholar

Lara, M. J., Genet, H., McGuire, A. D., Euskirchen, E. S., Zhang, Y., Brown, D. R. N., et al. (2016). Thermokarst rates intensify due to climate change and forest fragmentation in an Alaskan boreal forest lowland. Glob. Chang. Biol. 22, 816–829. doi: 10.1111/gcb.13124

PubMed Abstract | CrossRef Full Text | Google Scholar

Linell, K. A. (1973). “Long-term effects of vegetative cover on permafrost stability in an area of discontinuous permafrost. National Academy of Sciences, National Research Council.” in Proceedings of Permafrost: North American contribution to the Second International Conference. 13–28 July, 1973; 688–693.

Google Scholar

Luláková, P., Perez-Mon, C., Šantrůčková, H., Ruethi, J., and Frey, B. (2019). High-alpine permafrost and active-layer soil microbiomes differ in their response to elevated temperatures. Front. Microbiol. 10:668. doi: 10.3389/fmicb.2019.00668

PubMed Abstract | CrossRef Full Text | Google Scholar

Mackelprang, R., Waldrop, M. P., DeAngelis, K. M., David, M. M., Chavarria, K. L., Blazewicz, S. J., et al. (2011). Metagenomic analysis of a permafrost microbial community reveals a rapid response to thaw. Nature 480, 368–371. doi: 10.1038/nature10576

PubMed Abstract | CrossRef Full Text | Google Scholar

Masse, J., Prescott, C. E., Renaut, S., Terrat, Y., and Grayston, S. J. (2017). Plant community and nitrogen deposition as drivers of alpha and Beta diversities of prokaryotes in reconstructed oil sand soils and natural boreal Forest soils. Appl. Environ. Microbiol. 83, e03319–e03416. doi: 10.1128/AEM.03319-16

PubMed Abstract | CrossRef Full Text | Google Scholar

McMurdie, P. J., and Holmes, S. (2013). Phyloseq: An R package for reproducible interactive analysis and graphics of microbiome census data. PLoS One 8:e61217. doi: 10.1371/journal.pone.0061217

PubMed Abstract | CrossRef Full Text | Google Scholar

Messan, K. S., Jones, R. M., Doherty, S. J., Foley, K., Douglas, T. A., and Barbato, R. A. (2020). The role of changing temperature in microbial metabolic processes during permafrost thaw. PLoS One 15:e0232169. doi: 10.1371/journal.pone.0232169

PubMed Abstract | CrossRef Full Text | Google Scholar

Mickan, B. S., Abbott, L. K., Solaiman, Z. M., Mathes, F., Siddique, K. H. M., and Jenkins, S. N. (2019). Soil disturbance and water stress interact to influence arbuscular mycorrhizal fungi, rhizosphere bacteria and potential for N and C cycling in an agricultural soil. Biol. Fertil. Soils 55, 53–66. doi: 10.1007/s00374-018-1328-z

CrossRef Full Text | Google Scholar

Monteux, S., Weedon, J. T., Blume-Werry, G., Gavazov, K., Jassey, V. E. J., Johansson, M., et al. (2018). Long-term in situ permafrost thaw effects on bacterial communities and potential aerobic respiration. ISME J. 12, 2129–2141. doi: 10.1038/s41396-018-0176-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Neary, D. G., Klopatek, C. C., DeBano, L. F., and Ffolliott, P. F. (1999). Fire effects on belowground sustainability: a review and synthesis. For. Ecol. Manag. 122, 51–71. doi: 10.1016/S0378-1127(99)00032-8

CrossRef Full Text | Google Scholar

Numan, M. T., and Bhosle, N. B. (2006). α-L-Arabinofuranosidases: the potential applications in biotechnology. J. Ind. Microbiol. Biotechnol. 33, 247–260. doi: 10.1007/s10295-005-0072-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Oksanen, J., Blanchet, F. G., Friendly, M., Kindt, R., Legendre, P., McGlinn, D., et al. (2020). Vegan: Community Ecology Package. Available at: https://CRAN.R-project.org/package=vegan (Accessed March 16, 2022).

Google Scholar

Osburn, E. D., McBride, S. G., Aylward, F. O., Badgley, B. D., Strahm, B. D., Knoepp, J. D., et al. (2019). Soil bacterial and fungal communities exhibit distinct long-term responses to disturbance in temperate forests. Front. Microbiol. 10:2872. doi: 10.3389/fmicb.2019.02872

PubMed Abstract | CrossRef Full Text | Google Scholar

Parada, A. E., Needham, D. M., and Fuhrman, J. A. (2016). Every base matters: assessing small subunit rRNA primers for marine microbiomes with mock communities, time series and global field samples: primers for marine microbiome studies. Environ. Microbiol. 18, 1403–1414. doi: 10.1111/1462-2920.13023

PubMed Abstract | CrossRef Full Text | Google Scholar

Pardo, L. H., Hemond, H. F., Montoya, J. P., Fahey, T. J., and Siccama, T. G. (2002). Response of the natural abundance of 15N in forest soils and foliage to high nitrate loss following clear-cutting. Can. J. For. Res. 32, 1126–1136. doi: 10.1139/x02-041

CrossRef Full Text | Google Scholar

Partain, J. L., Alden, S., Strader, H., Bhatt, U. S., Bieniek, P. A., Brettschneider, B. R., et al. (2016). An assessment of the role of anthropogenic climate change in the Alaska fire season of 2015. Bull. Am. Meteorol. Soc. 97, S14–S18. doi: 10.1175/BAMS-D-16-0149.1

CrossRef Full Text | Google Scholar

Prosser, J. I., Head, I. M., and Stein, L. Y. (2014). “The Family Nitrosomonadaceae,” in The Prokaryotes. eds. E. Rosenberg, E. F. DeLong, S. Lory, E. Stackebrandt, and F. Thompson (Berlin, Heidelberg: Springer Berlin Heidelberg), 901–918.

Google Scholar

Quast, C., Pruesse, E., Yilmaz, P., Gerken, J., Schweer, T., Yarza, P., et al. (2012). The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 41, D590–D596. doi: 10.1093/nar/gks1219

PubMed Abstract | CrossRef Full Text | Google Scholar

R Core Team (2019). R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. Available at: https://www.R-project.org/ (Accessed March 16, 2022).

Google Scholar

Saidi-Mehrabad, A., Neuberger, P., Hajihosseini, M., Froese, D., and Lanoil, B. D. (2020). Permafrost microbial community structure changes across the Pleistocene-Holocene boundary. Front. Environ. Sci. 8:133. doi: 10.3389/fenvs.2020.00133

CrossRef Full Text | Google Scholar

Schütte, U. M. E., Henning, J. A., Ye, Y., Bowling, A., Ford, J., Genet, H., et al. (2019). Effect of permafrost thaw on plant and soil fungal community in a boreal forest: does fungal community change mediate plant productivity response? J. Ecol. 107, 1737–1752. doi: 10.1111/1365-2745.13139

CrossRef Full Text | Google Scholar

Schuur, E. A. G., Crummer, K. G., Vogel, J. G., and Mack, M. C. (2007). Plant species composition and productivity following permafrost thaw and Thermokarst in Alaskan tundra. Ecosystems 10, 280–292. doi: 10.1007/s10021-007-9024-0

CrossRef Full Text | Google Scholar

Schuur, E. A. G., McGuire, A. D., Schädel, C., Grosse, G., Harden, J. W., Hayes, D. J., et al. (2015). Climate change and the permafrost carbon feedback. Nature 520, 171–179. doi: 10.1038/nature14338

CrossRef Full Text | Google Scholar

Seitz, T. J., Schütte, U. M. E., and Drown, D. M. (2021). Soil disturbance affects plant productivity via soil microbial community shifts. Front. Microbiol. 12:619711. doi: 10.3389/fmicb.2021.619711

PubMed Abstract | CrossRef Full Text | Google Scholar

Sewell, P. D., Quideau, S. A., Dyck, M., and Macdonald, E. (2020). Long-term effects of harvest on boreal forest soils in relation to a remote sensing-based soil moisture index. For. Ecol. Manag. 462:117986. doi: 10.1016/j.foreco.2020.117986

CrossRef Full Text | Google Scholar

Shen, J., Chen, C. R., and Lewis, T. (2016). Long term repeated fire disturbance alters soil bacterial diversity but not the abundance in an Australian wet sclerophyll forest. Sci. Rep. 6:19639. doi: 10.1038/srep19639

PubMed Abstract | CrossRef Full Text | Google Scholar

Violle, C., Pu, Z., and Jiang, L. (2010). Experimental demonstration of the importance of competition under disturbance. Proc. Natl. Acad. Sci. 107, 12925–12929. doi: 10.1073/pnas.1000699107

PubMed Abstract | CrossRef Full Text | Google Scholar

Vitousek, P. M., Gosz, J. R., Grier, C. C., Melillo, J. M., Reiners, W. A., and Todd, R. L. (1979). Nitrate losses from disturbed ecosystems. Science 204, 469–474. doi: 10.1126/science.204.4392.469

PubMed Abstract | CrossRef Full Text | Google Scholar

Warnes, G. R., Bolker, B., Bonebakker, L., Gentleman, R., Liaw, W. H. A., Lumley, T., et al. (2020). Gplots: Various R Programming Tools for Plotting Data. Available at: http://CRAN.R-project.org/package=gplots (Accessed March 16, 2022).

Google Scholar

Wickham, H. (2016). ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag, New York. ISBN 978-3-319-24277-4. Available at: https://ggplot2.tidyverse.org (Accessed March 16, 2022).

Google Scholar

Wolken, J. M., Hollingsworth, T. N., Rupp, T. S., Chapin, F. S., Trainor, S. F., Barrett, T. M., et al. (2011). Evidence and implications of recent and projected climate change in Alaska’s forest ecosystems. Ecosphere 2:art124. doi: 10.1890/ES11-00288.1

CrossRef Full Text | Google Scholar

Wu, J., Wang, Y., Park, S.-Y., Kim, S. G., Yoo, J. S., Park, S., et al. (2016). Secreted alpha-N-Arabinofuranosidase B protein is required for the full virulence of Magnaporthe oryzae and triggers host defences. PLoS One 11:e0165149. doi: 10.1371/journal.pone.0165149

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, Y., Wu, L., Lin, Q., Yuan, M., Xu, D., Yu, H., et al. (2013). Responses of the functional structure of soil microbial community to livestock grazing in the Tibetan alpine grassland. Glob. Chang. Biol. 19, 637–648. doi: 10.1111/gcb.12065

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, Z., Wang, C., and Luo, Y. (2018). Effects of forest degradation on microbial communities and soil carbon cycling: a global meta-analysis. Glob. Ecol. Biogeogr. 27, 110–124. doi: 10.1111/geb.12663

CrossRef Full Text | Google Scholar