Biological invasions of plants have profound effects on ecosystem functioning by directly and indirectly altering soil microbiota, especially when invasive plants co-invade with their associated microbiomes. Ecosystem functions may recover slowly following invader removal, with implications for restoration.

We investigated the recovery of soil ecosystem function (measured as soil enzymes) following the removal, at different densities and times, of invasive Pinus spp. in New Zealand, and how different enzymatic activities responded to pine legacies.

Enzymatic activities were driven by pine legacies via both abiotic (soil nutrients) and biotic (fungi and bacteria) soil properties, with different enzymes showing distinct patterns. The activity of the enzymes cellobiohydrolase (cellulose degrading), β-glucosidase (cellulose degrading), N-acetyl-glucosaminidase (chitin degrading), laccase (lignin oxidising) and acid phosphatase (organic phosphate hydrolysing) were influenced by time since pine removal and by pine density at removal via effects on biotic communities. In comparison, Mn-peroxidase (lignin oxidising) was positively correlated with density of pines at removal and was negatively correlated with time since removal and was only influenced by fungal communities.

Synthesis. The recovery of soil enzymatic function following invasive species removal is slow, and dependent on pine legacies through the gradual changes in fungal and bacterial communities. The cascading effects of these changes suggest potential implications for the success of future plant establishment and restoration of co-invaded ecosystems.

At each site (22 sites in total), 5 soil samples were collected. From each sample, DNA was extracted for amplicon sequencing of fungal and bacterial communities. Soil nutrient concentrations and enzyme activities were also measured. Specific methodology is described below and in the manuscript.

DNA extraction

DNA was extracted from all soil samples using a DNeasyPowerSoil Pro Kit (Qiagen, Hilden, Germany) as per the manufacturer’s protocol. The fungal amplicon library was prepared in a one-step PCR using fungal specific primers fITS7 (5’-GTG ART CAT CGA ATC TTT G -3’) and ITS4 (5’-TCC TCC GCT TAT TGA TAT GC-3’) which amplify the ITS2 region (universal genetic barcode for fungi; Ihrmark et al. 2012). The ITS4 reverse primer included both Illumina adapted and index sequences (for identification of sequenced amplicons) and the fITS7 primer included Illumina adaptor sequences. The bacterial amplicon library was prepared in a one-step PCR using bacterial specific primers 515F (5’- GTG YCA GCM GCC GCG GTA A-3’) and 806R (5’-GGA CTA CNV GGG TWT CTA AT-3’) which amplify the 16S V4 region (Thompson et al. 2017). The 806R reverse primer included both Illumina adapted and index sequences (for identification of sequenced amplicons) and the 515F included Illumina adaptor sequences.

PCR products were visualised on 1% agarose gels and concentrations normalised using SequelPrep (Thermofisher). Uniquely indexed libraries (n = 110 for fungi, n = 110 for bacteria), along with positive and negative controls, were pooled for each sequencing run which were performed on an Illumina MiSeq using 2 x 250 paired end V2 chemistry following the manufacturer’s recommendations.

Bioinformatics analysis

Paired-end reads were merged using USEARCH v11 (Edgar 2010) with a minimum overlap length of 50 bp with no gaps allowed in the merged alignments. Sequence deconvolution such as quality control and clustering was carried out using VSEARCH (https://github.com/torognes/vsearch.wiki.git). Specifically, sequences less than 200 bp and low mean quality (<20) were removed. Sequences that passed quality control were clustered into operational taxonomic units (OTUs). Fungal OTUs were identified by searching against the public database provided by UNITE (UNITE Community 2017, Nilsson et al. 2018). Bacteria OTUs were identified by searching against the public database provided by SILVA (Quast et al. 2012).

Ecological information of each fungal OTU was based primarily on FUNGuild (https://github.com/UMNFuN/FUNGuild), following the user’s manual where assignments were made on ranks of ‘probable’ and ‘highly probable’ (Nguyen et al. 2016).

Abiotic soil properties

Soils were air dried and sieved at 2 mm in preparation for analysis. All tests were conducted at R J Hill Laboratories Limited (Hamilton, 3240, New Zealand). 

Soil enzyme assay

The enzyme activity of six enzymes were measured: cellobiohydrolase (cellulose degrading), β-glucosidase (cellulose degrading), N-acetyl-glucosaminidase (NAG; chitin degrading), acid phosphatase (organic phosphate hydrolysing), laccase (lignin oxidising) and Mn-peroxidase (lignin oxidising). These enzymes were chosen to encompass carbon, nitrogen and phosphorus cycling based on recent literature (German et al. 2011, Kyaschenko et al. 2017a, Sterkenburg et al. 2018). Cellobiohydrolase was assessed via the hydrolysis of 4-methylumbelliferyl β-D-cellobioside. β-glucosidase was assessed via the hydrolysis of 4-methylumbelliferyl β-D-glucopyranoside. NAG was assessed via the hydrolysis of 4-methylumbelliferyl N-acetyl-β-D-glucosaminide. Acid phosphatase was assessed via the hydrolysis of 4-methylumbelliferyl-phosphate. Soil samples were added to sodium acetate buffer (50 mM, pH 5) in the proportion of 1 g soil to 125 mL of buffer and the suspension was mixed for 1 min. For each enzyme, a 200 μL sample aliquot (pipetted during continuous stirring) was added to 50 μL of fluorogenic substrate. After incubation, 10 μL 1 M NaOH was added to stop the reaction and, after an additional 10 min, fluorescence was measured (excitation wavelength 365 nm, excitation slit 2.5, emission wavelength 450 nm and emission slit 2.5) and related to 4-methylumbelliferone standards.

Laccase activity was assessed by adding 100 μL of extract to 100 μL of reaction solution containing 50 mM 2,2'-azino-bis(3-ethylbenzthiazoline-6-sulphonic acid) (ABTS) in sodium acetate buffer. Peroxidase activity was assessed via the oxidative coupling of 3-(dimethylamino)benzoic acid (DMAB) and 3-methyl-2-benzothiazolinone hydrazine (MBTH) in the presence of Mn2+ and H2O2 (following protocols of Kyaschenko et al. 2017a). A 50 μL aliquot of extract was added to 150 μL of reaction solution (50 mM DMAB, 1 mM MBTH, 1 mM MnSO4, 100 mM sodium lactate, 100 mM sodium succinate, 67 μL H2O2, pH 5). Enzyme activities were estimated by regression of absorbance (laccase at 420 nm and Mn-peroxidase at 590 nm) against time. Boiled aliquots served as negative controls for all enzyme measurements.

All analyses were conducted in R version 4.1.0 (R Core Team 2021) as described in the manuscript.