DOI: 10.5061/dryad.ffbg79d33
title: Influence of soil type and vertical zonation on soil fungal communities associated with natural jack pine forests

description: Natural jack pine forests established on sandy eskers and clayey soils
 represent a unique habitat in the boreal forest in northeastern Canada and
 are generally associated with well-differentiated soil horizons, but
 little is known regarding their soil fungal diversity and how soil type
 and vertical zonation would shape their communities. To address this
 question, we characterized soil fungal communities from litter, organic
 and mineral horizons in 18 natural jack pine forests in esker and clayey
 soils using ITS2 rDNA metabarcoding. Ectomycorrhizal fungi dominated the
 organic and mineral horizons in esker soils, while saprotrophic fungi were
 dominant in clayey soils. The dominance of ectomycorrhizal fungi coincided
 with higher C accumulation within the organic horizon in esker soils. The
 nutrient-poor and more acidic esker sites harbored less diverse
 ectomycorrhizal and saprotrophic communities than clayey soils. Several
 stress-tolerant and melanin-rich fungal taxa were indicators of esker
 soils. However, fungal communities were more distinct between horizons
 than soil types and soil type effects differed along the soil profile. Our
 findings highlight the need to fill knowledge gaps in soil fungal
 diversity associated with jack pine forests and pave the way for adapted
 and multi-resource forest management in esker forest ecosystems
 considering their underground fungal specificities.
description: <strong>Study area</strong> The
 study was conducted in the Abitibi-Témiscamingue region, Quebec (Canada).
 The study area is located in northern Quebec’s clay plain, within the
 balsam fir-paper birch (<em>Abies balsamea</em> (Linnaeus)
 Miller<em>-Betula papyrifera</em>Marshall) bioclimatic domain
 (Grondin 1996). Climatic conditions are subpolar, subhumid, and
 continental with an annual average temperature and precipitation of
 respectively between 1-2 °C and 825-975 mm according to the Canadian
 climate normal 1981-2010 station data (Environment Canada, 2023). Forest
 vegetation mainly consisted of jack pine stands with interspersed stems of
 <em>Picea mariana</em> (Miller) Britton, Sterns &
 Poggenburgh, <em>Populus tremuloides</em> Michaux,
 <em>Betula papyrifera</em> Marshall, <em>Abies
 balsamea</em> (Linnaeus) Miller, and <em>Acer
 rubrum</em> Linnaeus. Specifically, the Abitibi region is
 characterized by natural gradients of soil physicochemical properties,
 from sandy eskers to clayey soil types. Indeed, esker soils account for
 ~3% of the Abitibi’s surface while the remainder of the region’s flat
 topography (~57%) is dominated by clayey soils (Robitaille and Saucier
 1998; Cloutier et al. 2007). Abitibi eskers have mineral surface deposits
 with an overall thickness of more than 25 cm, an illuvial horizon with
 coarse texture, and a high stone content. Conversely, glacio-lacustrine
 clayey soils, which are associated with Luvisols or Gleysols in the
 Abitibi region (Laverdière and De Kimpe 1984; Bergeron et al. 2007), are
 characterized by an illuvial horizon with medium texture and a low stone
 content (Blouin and Berger 2002). <strong>Site
 selection</strong> Site selection was based on
 four main criteria, namely soil type (sandy esker/clayey), dominant
 vegetation (jack pine, minimum 50-75% of canopy cover), stand age (minimum
 50-60 years), stand origin (natural fire-origin, totally or mostly
 unimpacted by anthropogenic disturbance), and accessibility (maximum 350 m
 from forest roads). We used satellite images and forest inventory data
 from the Quebec government’s online platform ‘Forêt Ouverte’ for site
 selection. We selected 18 natural jack pine stands (12 on sandy eskers, 6
 on the clay plain). Only six clayey sites fitting our selection criteria
 were found within the prospected territory due to the scarcity of
 unmanaged (natural) jack pine stands on clayey soils. Sites were separated
 by a minimum distance of 5 km to ensure independence between samples and
 avoid pseudoreplication. Finally, the study design consisted of a
 first-level factor (soil type) with two levels (sand and clay), and a
 second-level factor (soil horizons) with three levels (litter, organic,
 and mineral). <strong>Soil
 sampling</strong> Soil samples were collected
 during the summer of 2022 to characterize soil fungal communities and
 obtain the physicochemical composition across the soil profile. A circular
 plot with a 28 m radius and an area of 2,500 m<sup>2</sup> was
 established at each study site, at least 20 m from the forest edge to
 avoid related variability bias (Dickie and Reich 2005). Within each plot,
 five jack pine trees were randomly selected (at least 10 m apart), and two
 soil sample replicates per tree (2-3 m apart) were collected, following a
 modified protocol by Tedersoo et al. (2014; 2021). Sampling locations
 around trees were randomly selected with the restriction criteria of being
 strictly opposite (angle of 180°). In total, 180 soil cores (2 replicates
 × 5 trees × 18 sites) were collected. For each sampling
 spot, forest fallen litter (dead leaves, needles) (hereafter referred to
 as litter horizon) was collected in plastic bags for DNA and
 physicochemical analyses, followed by samples from the organic layer
 (hereafter referred to as organic horizon). The latter represented a
 dark-colored layer rich in organic matter at various decomposition stages,
 mainly composed of fallen plant material. These organic samples included
 both F and H horizons. After the removal of the organic layer, we sampled
 the mineral soil (hereafter referred to as mineral horizon) using a
 pedological auger (25 cm in length and 7.5 cm in diameter), capturing both
 eluvial and illuvial horizons. For esker sites, the mineral soil profile
 was characterized by the formation and accumulation of organo-metallic
 assemblages, involving Al and Fe elements, in a leached grey eluvial and
 rust-brown illuvial horizon, respectively. Due to the variable thickness
 of the eluvial horizons within both esker and clayey soil profiles, the
 proportion of this horizon in the mineral soil samples varied among sites.
 Large roots and coarse woody debris were systematically removed from
 organic and mineral material while sampling. Samples
 were pooled per horizon for each site, resulting in one composite sample
 for each horizon. This represents 18 composite samples per horizon for a
 total of 54 composite soil samples (3 horizons × 18 sites). Each composite
 soil samples were divided into two sub-samples: one for eDNA-based soil
 fungal community analysis and the other for physicochemical analysis.
 Composite organic and mineral soil samples for eDNA analysis were sieved
 with a 6 mm mesh in the field, placed in plastic bags, transported in an
 ice-filled cooler and stored at – 25°C in the Ecology Research Group of
 Abitibi RCM (GREMA) laboratory until further processing. Composite litter
 samples for eDNA analysis were crushed in liquid nitrogen using a mortar
 and pestle prior to DNA extraction. Composite soil samples for
 physicochemical analyses were sieved with a 4 mm mesh using an automatic
 vibrating sieve AS 200 Control (ATS Care Retsch, Haan, Germany) after
 being forced-air dried at room temperature for 14 days to facilitate the
 sieving process. <strong>Soil physicochemical
 analyses</strong> Composite soil samples
 (&gt; 50g) were sent to the organic and inorganic chemistry laboratory
 of the Forest Research Direction (Quebec, QC, Canada) for physicochemical
 analysis. Total nitrogen (N) and carbon (C) contents were determined by
 combustion using a CN 928 elemental analyzer (LECO Corp., St Joseph, MI,
 USA) with thermal conductivity detection for nitrogen and non-dispersive
 infra-red (NDIR) cell detection for carbon. Organic matter content
 (hereafter OM) and percentage of humidity (hereafter humidity) were
 determined by incineration, a method commonly referred to as the loss on
 ignition (LOI) method (Davies 1974). Soil pH was measured using 10 g of
 soil mixed with 20 mL of distilled water with an Orion VersaStar Pro pH
 meter (Thermo Fisher Scientific Inc., Pittsburgh, PA, USA). Elements (P,
 K, Ca, Mg, Mn, Zn, Al, Fe) were extracted using the Mehlich-III method
 (Mehlich 1984) and measured by plasma atomic emission spectroscopy (Optima
 8300 model, ICP-OES, Perkin Elmer, Waltham, MA, USA). Particle size
 distribution and textural class determination were performed on soil
 samples consisting of fine earth (&lt;2 mm) containing 5% or less
 carbon using the Bouyoucos method (Bouyoucos 1962).
 <strong>DNA extraction, amplification and library
 preparation</strong> DNA extraction was carried
 out at the Environmental Genomic Laboratory of the Laurentian Forestry
 Centre (EGL-CFL) on 150 mg of composite litter and organic samples, and
 250 mg of composite mineral samples using a DNeasy powersoil pro kit
 (Qiagen, Valencia, CA, US) and the QIAcube automated instrument (Qiagen,
 Hilden, Germany) in accordance with the manufacturer’s instructions. DNA
 was quantified with the Qubit™ dsDNA BR Assay Kit (Thermo Fisher
 Scientific Inc, Wilmington, USA). Fungal DNA amplicon libraries were
 prepared at the EGL-CFL following a previously described procedure (Samad
 et al. 2023) and sequenced on an Illumina MiSeq platform at the Next
 Generation Sequencing Platform of the CHU de Québec-Université Laval
 Research Centre using a MiSeq v3 600-cycle Reagent Kit. Even though one
 ITS region contains less genetic information than the entire ITS (Tedersoo
 et al. 2022), the ITS2 region of the fungal ribosomal DNA was amplified
 using the primer set ITS9F (5’-GAACGCAGCRAAIIGYGA-3’) and ITS4R
 (5’-TCCTCCGCTTATTGATATGC-3’ (White et al. 1990; Ihrmark et al. 2012;
 Rivers 2016) because of limitations to sequence length that can be
 obtained with Illumina Sequencing. Negative controls in the DNA extraction
 and PCR amplification were used to control potential contaminants during
 the process. <strong>Bioinformatic and sequence
 analyses</strong> All bioinformatics analyses
 were performed in QIIME2 v2023.2.0 (Bolyen et al. 2019). Demultiplexed
 FASTQ (R1 and R2) files were first imported in QIIME2 using the ‘qiime
 import’ command which converts the data into a QZA archive file so it can
 be processed by the rest of the QIIME2 workflow. First, primers were
 removed using the QIIME2 implementation of CutAdapt (Martin 2011).
 Resulting forward and reverse reads went through the DADA2 (Callahan et
 al. 2016) pipeline for sequence quality control and feature table
 construction using the ‘’qiime dada2 denoise-paired’’ command. During this
 process, low-quality regions of the sequences were trimmed, paired reads
 assembled, chimeric sequences filtered, remaining high-quality sequences
 dereplicated, and singletons and very low-frequency abundance ASVs removed
 using the “qiime feature-table filter-features” command. The output of
 this step was a feature table (i.e., ASV table) which contains read counts
 for each unique sequence in each sample of the dataset, and feature data
 which contains sequences corresponding to each ASV. A
 taxonomy was assigned to each ASV based on the UNITE reference database
 (version 9.0) (Abarenkov et al. 2010, 2023; Kõljalg et al. 2005, 2013;
 Nilsson et al. 2019b) using the ‘‘classify-sklearn’’ with a Naive Bayes
 classifier. All ASVs not assigned to the ‘Fungi’ kingdom were removed from
 further analyses using the “qiime taxa filter-table”. The resulting fungal
 ASV table was used to perform downstream statistical analyses. Functional
 annotation of ASVs was performed for most fungal guilds using FUNGuild
 version 1.1 (Nguyen et al. 2016).
description: # Influence of soil type and vertical zonation on soil fungal communities
 associated with natural jack pine forests
 [https://doi.org/10.5061/dryad.ffbg79d33](https://doi.org/10.5061/dryad.ffbg79d33) Three different datasets are presented. **First**, an XLSX file (**Fungal_ASV_dataset_non_rarefied.xlsx**) containing the occurrence of fungal ASVs across our samples, including extraction and PCR controls, based on an ITS2 Illumina metabarcoding (environmental DNA) technique. This dataset is non rarefied. This dataset also includes the taxonomic classification for each fungal ASV (at the Class, Order, Family, Genus, and Species levels), as well as ecological classification and attributes (Trophic Mode, Guild, Growth Morphology, Trait, etc.). The confidence of these assignations is also indicated, as well as the source if applicable. **Second**, an XLSX file (**Fungal_sequence_correspondence_dataset.xlsx**) containing the correspondence between the ASVs presented in the first dataset (column ASV_ID) and the associated sequences (column Sequence). **Third**, an XLSX file (**Chemistry_dataset.xlsx**) containing chemical properties for our samples. The enumeration of the samples mirrors the one from the first dataset. **Descriptions** **Fungal_ASV_dataset_non_rarefied.xlsx** * *ASV_ID*: Name (or identity) of the different ASVs. * *ITS2-CTRL-neg-extraction0*: Negative control sample for DNA extraction, number 0. * *ITS2-CTRL-neg-extraction1*: Negative control sample for DNA extraction, number 1. * *ITS2-CTRL-neg-extraction2*: Negative control sample for DNA extraction, number 2. * *ITS2-CTRL-neg-PCR*: Negative control sample for PCR amplification. Sample indexation follows this structure: ITS2-(Soil type)(Site number)(Horizon). For soil type, it can be either SE = esker soil, or SA = clayey soil. For esker soil, site number ranges from 1 to 12, and 1 to 6 for clayey soil. For horizons, it can be either L = litter horizon, O = organic layer, or M = mineral soil. The combination of these three parameters made up a unique index for each sample of our study. For exemple, sample *ITS2-SA1O* corresponds to the organic layer from clayey soil of site 1. Sample are listed from column 6 to 59 (total of 54 samples). The taxonomic classification fo each ASV, based on the FUNGuild database, is listed according to seven levels: *Kingdom, Phylum, Class, Order, Family, Genus*, and *Species*. The taxonomic classification follow this structure: (First letter of the classification level)(Taxonomy). For the first letter of the classification level, it can be either k = Kingdom, c = Class, o = Order, f = Family, g = Genus, and s = Species. If we failed to assign a taxonomic classification at a given level, it is indicated as (First letter of the classification level)__unidentified. Since ASVs were filtered for Fungi, the column Kingdom contains only one category (K__Fungi). ASVs were also assigned an "ecological" classification, based on the FUNGuild database, listed as *Trophic Mode*, *Guild* and *Growth Morphology*. The detailed definitions indicated below refers to the information provided by FUNGuild ([https://github.com/UMNFuN/FUNGuild](https://github.com/UMNFuN/FUNGuild) * *Trophic Mode*: One of three trophic categories [Pathotroph = receiving nutrients at the expense of the host cells and causing disease (e.g., biotroph, parasite, pathogen, etc.); Saprotroph = receiving nutrients by breaking down dead host cells (e.g., wood rotters, litter rotters, etc.); Symbiotroph = receiving nutrients by exchanging resources with host cells (e.g., ectomycorrhiza, lichens, etc.)]. * *Guild*: Provide a relevant guild descriptor [Pathotroph: Animal Pathogen (including human pathogens - typically annoated as such), Bryophyte Parasite, Clavicipitaceous Endophyte, Fungal Parasite, Lichen Parasite, Plant Pathogen; Saprotroph: Dung Saprotroph (i.e., coprophilous), Leaf Saprotroph (e.g., leaf litter decomposer), Plant Saprotroph, Soil Saprotroph (e.g., rhizosphere saprobe - typically annoated as a rhizosphere fungus), Undefined Saprotroph (e.g., a general saprobe, or in cases where the ecology is not known but suspected to be a saprobe), Wood Saprotroph (e.g., wood rotting fungi); Symbiotroph: Ectomycorrhizal, Ericoid Mycorrhizal, Endophyte, Epiphyte, Lichenized (i.e., lichen)]. * *Growth Morphology*: Basic morphological categories including "Agaricoid", "Boletoid", "Clathroid", "Clavarioid", "Corticioid", "Facultative Yeast", "Gasteroid", "Hydnoid", "Microfungus", "Phalloid", "Polyporoid", "Rust", "Secotioid", "Smut", "Thallus", "Tremelloid", or "Yeast" should be noted. * *Trait*: Other functional or morphological traits such as "Brown Rot", "Hypogeous", "Poisonous", "Soft Rot", "White Rot" or "Contact Exploration Type" would be appropriate for this field. * *Confidence Ranking*: "Highly Probable" (= absolutely certain), "Probable" (= fairly certain), "Possible" (= suspected but not proven, conflicting reports given, etc.). * *Notes*: Any other relevant information related to the taxon (e.g., "Facultative human pathogen causing coccidioidomycosis" for Coccidioides immitis, "Host - Arecaceae - Casuarinaceae, Palmae, Sapotaceae" for Xylaria obovata, or "Syn. Umbilicaria" for Actinogyra). * *Citation/Source*: Publication, website, etc. from which the corresponding guild information was derived. Data for taxa that is based on peer-reviewed publications is preferred. **Fungal_sequence_correspondence_dataset.xlsx** * *ASV_ID*: Same ASV name (or identity) as in **Fungal_ASV_dataset_non_rarefied.xlsx** * *Sequence*: Generated ITS2 sequence for each ASV. **Chemistry_dataset.xlsx** * *Site_number*: Same site indexation as in **Fungal_ASV_dataset_non_rarefied.xlsx**. Site number for esker soils ranges from 1 to 12, and from 1 to 6 for clayey soils.  * *Soil_type*: Esker = esker soil, Clay = clayey soils. * *Horizon*: Litter = litter horizon, Organic = organic layer, M = mineral soil. * *Organic_matter*: Organic matter content (g/Kg) * *Humidity*: Percentage of humidity (%). * *Total_C*: Total carbon content (g/Kg) * *Total_N*: Total nitrogen content (g/Kg) * *P*: Phosphorus content (mg/Kg) * *K*: Potassium content (mg/Kg) * *Ca*: Calcium content (mg/Kg) * *Mg*: Magnesium content (mg/Kg) * *Mn*: Manganese content (mg/Kg) * *Zn*: Zinc content (mg/Kg) * *Al*: Aluminium content (mg/Kg) * *Fe*: Iron content (mg/Kg) * *pH*: Potential of hydrogen.