Step 1: Sequence processing
Before processing, sequences with low quality reads were removed. ASV sequence counts were converted to relative abundances to standardize across samples with different sequencing depths. After processing, 7279 unique ASVs were kept.
Table 1. Raw data table.
Each row in the dataset represents a single sample. Predictor variables include tree species (pine, spruce, aspen), burn severity (control, low, medium, high), and year (2023, 2024). Response variables include ASV sequence composition and fungal order composition. Samples were categorized into seven treatment groups based on year and severity: 2023 (control, low, medium, high) and 2024 (low, medium, high).
Table 2. ASV sequences relative abundance table.
Column names were ASV sequences and row names are sample ID. Cells are relative abundances. Bray–Curtis dissimilarity was calculated from sample-level ASV relative abundances, and community differences were visualized using PCoA. Relative abundances were not averaged within treatment.
Step 2: Rough community composition comparison
Fire effect on communities
Figure 2. (a) PCoA illustrating differences in sapwood fungal community composition of lodgepole pine across low, medium, and high burn severity classes following the 2023 wildfire. Ellipses represent 95% confidence intervals around group centroids. (b) Dispersion plot illustrating within-group variation.
Figure 3. (a) PCoA illustrating differences in sapwood fungal community composition of white spruce across low, medium, and high burn severity classes following the 2023 wildfire. Ellipses represent 95% confidence intervals around group centroids. (b) Dispersion plot illustrating within-group variation.
Figure 4. (a) PCoA illustrating differences in sapwood fungal community composition of trembling aspen across low, medium, and high burn severity classes following the 2023 wildfire. Ellipses represent 95% confidence intervals around group centroids. (b) Dispersion plot illustrating within-group variation.
perMANOVA summary statistics for fire effects
Lodgepole pine
· C2023-low2023 =insignificant
· C2023-med2023= insignificant
· C2023-high2023 = significant <0.001*** (PERMANOVA, pseudo-F1,23=2.44, p<0.001)
Spruce
· C2023-low2023 = insignificant
· C2023-med2023 = insignificant
· C2023-high2023=significant <0.001*** (PERMANOVA, pseudo-F1,27=1.89, p<0.001)
Aspen
· C2023-low2023 =significant <0.05* (PERMANOVA, pseudo-F1,23=1.44, p<0.05)
· C2023-med2023 =significant <0.01** (PERMANOVA, pseudo-F1,28=1.72, p<0.05)
· C2023-high2023=significant 0.05*(PERMANOVA, pseudo-F1,29=1.68, p<0.01)
Changes in fungal communities after one year
Figure 5. PCoA illustrating change in sapwood fungal community composition of lodgepole pine from 2023 to 2024 across severity classes. Ellipses represent 95% confidence intervals around group centroids.
Figure 6. PCoA illustrating change in sapwood fungal community composition of white spruce from 2023 to 2024 across severity classes. Ellipses represent 95% confidence intervals around group centroids.
Figure 7. PCoA illustrating change in sapwood fungal community composition of trembling aspen from 2023 to 2024 across severity classes. Ellipses represent 95% confidence intervals around group centroids.
perMANOVA summary statistics for changes over time
Lodgepole pine
· Low2023-low2024 = significant <0.001*** (PERMANOVA, pseudo-F1,27=1.54, p<0.001)
· Med2023-med2024= significant <0.001*** (PERMANOVA, pseudo-F1,28=2.23, p<0.001)
· High2023-high2024= significant >0.01** (PERMANOVA, pseudo-F1,23=1.97, p<0.01)
Spruce
· Low2023-low2024= significant <0.01**
· Med2023-med2024= significant <0.001*** (PERMANOVA, pseudo-F1,30=1.94, p<0.01)
· High2023-high2024 =insignificant
Aspen
· Low2023-low2024 =insignificant
· Med2023-med2024 =insignificant
· High2023-high2024 = significant <0.05* (PERMANOVA, pseudo-F1,28=1.49, p<0.05)
Step 3: community composition comparison at the fungal order level
Table 3. Taxonomic assignment based on ASV sequences.
After Bray-Curtis, individual ASV sequences were assigned a taxonomy, which included its fungal order assignment. Trophic modes (also known as feeding lifestyles) were also assigned. Table 2 and Table 3 were combined to get the relative abundances of each fungal order in each sample. Relative abundances were then averaged within treatment groups and compared via stacked bar charts (see results page).
Limitations
The dataset is high-dimensional and sparse, with many low abundance ASVs and zero values. Sample sizes are not perfectly balanced across treatment groups due to sampling error in the field, missing samples, or removed due to poor sequencing quality. These characteristics influence the dispersion and statistical interpretations in downstream analysis. Within-group dispersion is observed to vary greatly among treatments potentially affecting the interpretation of both p-values and R² (“effect size”). perMANOVA is best for handling multivariate datasets like this, especially with multiple levels within each variable (E.g., burn severities: low, medium, high).
Another limitation was that there were no controls collected in 2024. This was intentional, because collecting one year after fire event will not represent true undisturbed conditions. Instead, we compared 2024 samples to 2023 samples to evaluate whether one-year-old communities converge toward, or diverge from, pre-disturbance conditions.