Rapid Microbial Community Changes During Initial Stages of Pine Litter Decomposition

Plant litter decomposition is a process enabling biogeochemical cycles closing in ecosystems, and decomposition in forests constitutes the largest part of this process taking place in terrestrial biomes. Microbial communities during litter decomposition were studied mainly with low-throughput techniques not allowing detailed insight, particularly into coniferous litter, as it is more difficult to obtain high quality DNA required for analyses. Motivated by these problems, we analyzed archaeal, bacterial, and eukaryotic communities at three decomposition stages: fresh, 3- and 8-month-old litter by 16/18S rDNA pyrosequencing, aiming at detailed insight into early stages of pine litter decomposition. Archaea were absent from our libraries. Bacterial and eukaryotic diversity was greatest in 8-month-old litter and the same applied to bacterial and fungal rDNA content. Community structure was different at various stages of decomposition, and phyllospheric organisms (bacteria: Acetobacteraceae and Pseudomonadaceae members, fungi: Lophodermium, Phoma) were replaced by communities with metabolic capabilities adapted to the particular stage of decomposition. Sphingomonadaceae and Xanthomonadaceae and fungal genera Sistotrema, Ceuthospora, and Athelia were characteristic for 3-month-old samples, while 8-month-old ones were characterized by Bradyrhizobiaceae and nematodes (Plectus). We suggest that bacterial and eukaryotic decomposer communities change at different stages of pine litter decomposition in a way similar to that in broadleaf litter. Interactions between bacteria and eukaryotes appear to be one of the key drivers of microbial community structure. Electronic supplementary material The online version of this article (10.1007/s00248-018-1209-x) contains supplementary material, which is available to authorized users.


Introduction
Decomposition of organic matter is a key step in nutrients cycling in all ecosystems. Plant litter decomposition is an important ecological process enabling biogeochemical cycles closing in terrestrial ecosystems. Decomposition in forests constitutes the largest part of this process taking place in terrestrial biomes, due to their immense area (~30% of land surface [1]) and large quantities of organic matter stored.
Forest litter decomposition was extensively studied since the early 1970s (reviewed, e.g., in [2]). However, chemical changes in decomposing materials as well as element cycling in relation to temperature and precipitation were assessed, while microbial and macroorganismal aspects of this process were less intensively studied [3,4]. Nevertheless, it is estimated that microbes are responsible for up to 90% of organic matter decomposition [5], and the dominating primary decomposers in boreal and temperate forest soil systems are microorganisms, mainly fungi and bacteria. The structure and development of decomposer communities can influence the pattern of decay [3].
As the material decomposes, chemical composition of the litter changes, and there is a shift from carbohydrates and aliphatic components constituting the largest pool in initial litter, to aromatic compounds at late stages of decomposition.
Electronic supplementary material The online version of this article (https://doi.org/10.1007/s00248-018-1209-x) contains supplementary material, which is available to authorized users. DNA and spectrum of lysed organisms. We hypothesized that the treatment would increase the yield and diversity of 16/18S rDNA amplicons. As there were no reports of Archaea being found in litter, we expected that they would be absent from pine litter samples. Bacterial and eukaryotic communities would be different at different stages of decomposition, specifically (i) microbial diversity would follow diversity of substrates, i.e., would be lower at later stages of decomposition, (ii) late communities would be more dominated by specialists capable of lignocellulose and lignin utilization, (iii) phyllosphere-related organisms would prevail initially and then would be replaced by those coming from soil, and (iv) a shift in metabolic capabilities of the community was also expected wherein organisms utilizing soluble small molecules initially present in the phyllosphere would be superseded by those in whose genomes reside genes enabling cellulose and, later, lignin utilization. To test these hypotheses, we prepared, sequenced, and analyzed pyrosequencing (454) libraries of archaeal and bacterial 16S, as well as eukaryotic 18S rRNA gene fragments derived from litter samples at three stages of decomposition.

Study Design and Samples
The study site was located in Sierbowice, southern Poland (GPS coordinates: N 50°34′ 15.20″, E19°39′ 55.40″). The vegetation at the site is dominated by approx. 50 years old Scots pine (Pinus sylvestris L.) with a small admixture of birch (Betula pendula Roth) and artificially introduced red oak (Quercus rubra L.), with European blueberry (Vaccinium myrtillus L.) and mosses in the groundcover.
We were interested in investigating three time points of initial stages of litter decomposition: t 0 = no field incubation, t 1 = 92 days, and t 2 = 242 days of field incubation, because we expected changes in community composition in line with rapid initial mass loss. In order to do so, freshly fallen brown pine litter from over a dozen of randomly picked trees in the area of about 3 ha was collected in October 2012. The needles were obtained by shaking the selected trees and their branches. A PVC foil was spread under the trees to facilitate litter collection and to prevent contamination of fresh litter with microbial communities from soil.
The collected litter was transported to the laboratory in plastic bags. The material was thoroughly mixed and split into two major parts. The first part consisted of three subsamples dedicated for water content analysis and DNA extractions to determine the t 0 microbial community. The second major part of the collected material was air dried in room temperature for a month. Three subsamples from the air-dried material were taken for chemical analyses, and the rest was used to make litter bags. The litter bags were made out of a nylon mesh (20 × 20 cm, mesh size 1 × 1 mm) and were packed with2 0 g of air dry litter. On the 1 December 2012, they were placed back in the field under the soil organic layer. Six bags were collected at each of the time points. Every time, three litter bags were dedicated for chemical and the remaining ones for molecular analyses. Before the DNA extractions, the samples were placed in a climatic chamber and acclimated for a week at 22°C to 70% of water holding capacity (WHC) and frozen at − 80°C. The acclimation step was performed in order to standardize the physiological state of microorganisms in the litter, as the collection of litter took place during different seasons. In the end, we obtained three biological replicates per time point, and three technical replicates (independent isolations) were made for each of four isolation methods used. Just as in case of the t 0 series before the chemical analyses, t 1 and t 2 samples were air dried.

Physicochemical Analyses
Water content was measured gravimetrically in fresh litter samples, immediately after their transfer to the laboratory. pH was measured in a slurry of 0.5-g air dried litter in 15 ml of demineralized water. The concentrations of particular elements were measured in powdered and dried material (12 h, 105°C). C and N concentrations were determined with Vario EL Cube (Elementar, Germany). The remaining macroelements and microelements were determined with a PinAAcle 900 Z atomic absorption spectrometer (Perkin Elmer, USA) after wet digestion with nitric acid in Titan MPS microwave sample preparation system (Perkin Elmer, USA).

Climatic Data
Climatic data (daily averages) were downloaded from dane. imgw.pl. As there is no meteorological station in the immediate vicinity of the sampling area, values from six nearest stations (Silniczka, Lgota Górna, Katowice-Pyrzowice, Olewin, Miechów, Jędrzejów-Sudół) were averaged and plotted in R.

DNA Isolation
PowerSoil DNA Isolation kit (MoBio, USA) was used according to the producer's protocol for PowerLyzer 24 bead beater involving one 45-s cycle of bead beating at 4000 rpm for control isolations and with the following modifications: (i) addition of lysozyme (Sigma, USA) to the final concentration of 2 mg/ml and achromopeptidase (500 U/ml, Sigma, USA) to the C1 buffer and incubation at 37°C for 1 h (AL method), (ii) addition of chitinase (Sigma, 0.01 U/ml) and incubation at for 1 h (Ch method) and (iii) addition of lisozyme, achromopeptidase and chitinase, and incubation as above (ALCh method). Isolated DNA was quantified with Qubit HS DNA kit (Invitrogen, USA), and the quality was measured spectrofotometrically on NanoDrop ND-1000 (Thermo Fisher Scientific, USA). DNA content was expressed in ng per 1 g of fresh litter.

Primers Design
Primers were designed basing on SILVA v.119 alignment [36]. The alignment was split into kingdom-specific parts with Mothur's get.lineage, and consensus sequences were generated at the 97% identity level (consensus.seqs). Visual inspection of resulting summary files allowed identification of highly conserved regions. Candidate pairs were checked with the online TestPrime tool [37] and IDT Oligoanalyzer [38]. Pinespecific primers were designed in the same way, but an alignment of Pinus 18S rRNA sequences from SILVA was used.

16S and 18S rRNA Gene Fragments Amplification and Pyrosequencing
Libraries of 16S/18S rRNA gene fragments were created with the use of two-step method, involving gene-specific primers tagged with M13/M13R overhangs in the first round of PCR and M13 bearing 9-nt MID sequences [39] and A adapter sequence (Roche, Switzerland) paired with M13R with B adapter overhang in the second round. Primer sequences and PCR conditions are listed in Table 1. The final products were quantified with Qubit HS DNA kit (Invitrogen), and 36 of them were pooled in equimolar amounts for each library. The library quality was assessed with HS DNA chip on Bioanalyzer (Agilent, USA). Libraries were emPCR amplified with the use of Titanium Lib-L kits (Roche, Switzerland) and sequenced on GS-Junior machine with Titanium chemistry (Roche, Switzerland) as per the manufacturer's protocols.

qPCR
Real-time PCR analyses were conducted using primers listed in Table 1 and FastStart SYBRGreen kit (Roche, Switzerland) on LightCycler 480 machine (Roche, Switzerland). The reaction mix included 3 pmol of forward and reverse primers, 2 ng of template DNA, 5 μl of 2 x concentrated kit, and water up to 10 μl. Standards were prepared from pure amplicons generated with primer pairs used for qPCR on DNA isolated from Escherichia coli, P. sylvestris, and Boletus badius. Standard curves were replicated five times, and samples were assayed in triplicates. Each run included negative control (water). Resulting numbers of copies were converted to copies/g of litter.

Bioinformatics Analyses
Pyrosequencing reads were processed with Mothur v. 1.32 [40] and custom-tailored Perl scripts, with modifications increasing the aggressiveness of denoising, chimera removal, as well as producing ten subsamples of the whole data and averaging the shared OTU table over those subsamples, as described earlier [41,42]. Brief summary of the procedure is given below. The flows were extracted from the .sff files, forward, and reverse reads separately (sffinfo), then they were assigned to samples basing on the MID sequences, trimmed to min. 500 and max. 650 flows (trim.flows), and denoised with AmpliconNoise algorithms (shhh.flows and shhh.seqs; [43]). Primers and MIDs were removed from the denoised seuqences, and the read set was dereplicated (unique.seqs) and aligned to the SILVA v.119 template alignment (align.seqs). Reads covering the desired region of the alignment (pos. 6500-22,500 for bacteria and 13,876-22,550 for eukaryotes) were chosen (screen.seqs) and gap only, and terminal gap-containing columns were removed from the alignment (filter.seqs). The set was dereplicated again, and residual sequencing and PCR noise was removed with Single Linkage pre-clustering (pre.cluster; [44]). Chimera identification and removal were performed in two rounds: (i) with UCHIME (chimera.uchime; [45]) and (ii) with PERSEUS (chimera.perseus; [43]).
Full-length sequences (list.seqs, get.seqs) were used for classification with naive Bayesian classifier (classify.seqs; [46]) using SILVA 119 template and taxonomy files (http:// www.mothur.org/w/images/2/27/Silva.nr_v119.tgz, accessed on September 4, 2014) for classification of bacterial reads and PR2 database [47] for eukaryotic ones at the bootstrap confidence level of 80%. Taxons Bunknown^and, in case of bacterial data, Bchloroplast^were removed from the final set. OTUs at the 0.03 dissimilarity level were constructed via average linkage (UPGMA), and singletons together with doubletons were removed from the data (remove.rare).
For the initial analyses of enzyme influence, the final read set was subsampled to 500 reads per sample ten times (sub.sample and regular expressions in the sed editor), subsamples were combined (cat), the whole set was dereplicated and used for distance matrix calculation (dist.seqs), and OTU construction via average neighbor clustering at 97% similarity level (cluster). Shared OTU table was constructed (make.shared), and averaged table was calculated with a Perl script (average_shared.perl). Final analyses were performed on a dataset in which all enzymatic treatments coming from one sample were combined. This dataset was subsampled to 700 (bacteria) and 1000 (eukaryota) reads per sample and processed as described above.
Relaxed neighbor joining (RNJ) tree was constructed from the final alignment with clearcut (clearcut; [48]). UniFrac [49] distance matrices were calculated in Mothur (unifrac.unweighted and unifrac.weighted) with subsampling the RNJ tree to include 700 and 1000 reads per sample for bacteria and eukaryota, respectively. Morisita-Horn [50] and Bray-Curtis [51] distance matrices were calculated in R (vegdist). NMDS and CCA analyses were performed in R with vegan's [52] metaMDS and cca functions, respectively. For NMDS, 1000 tries were used, and the same number of permutations was adopted in CCA. CCA models were built by backward selection with ordistep.
Co-correspondence analysis [53] was performed with the use of the coca function from the cocorresp R package [54]. Significance of the extracted axes was tested with permutation test (permutest), while the percent of fit was checked by leaveone-out cross-validation (crossval).

PICRUSt Analysis
For PICRUSt [55] analysis, bacterial sequences were classified as described above, but using GreenGenes taxonomy files (v. 13_8; [56]). Taxonomic information along with OTU table was converted to a .biom file using Mothur's make.biom function. The file was then converted to v.1.0.0 format using biom convert, as per https://github.com/rprops/PICRUSt_from_ mothur (visited Feb. 15, 2017). Normalized OTU table was generated with normalize_by_copy_number.py, and predicted metagenomes as well as NSTI scores were calculated with predict_metagenomes.py. Functions were collapsed to pathways at level 2 using categorize_by_function.py. Pathways pertaining to organismal systems were considered spurious and removed. Predicted metagenomes were analyzed with STAMP [57].

Statistical Analyses
R [58] was used for statistical computations, with Hmisc [59], phyloseq [60] and vegan [52] packages. Significance level used was 0.01. ANOVA with Tukey's HSD test was used to test for significance of differences in means, PERMANOVA (vegan's adonis), AMOVA (amova of ade4 package), as well as ANOSIM (vegan's anosim) for testing separation of clusters. Significance level used was 0.01. When testing for significance of grouping by enzymatic treatments, permutations were restricted to a given decomposition level (strata = decomposition). PERMDISP test (vegan's betadisper) was used to test the homogeneity of variance in community data.
Differences in KEGG pathway content between litter decomposition stages were tested in STAMP using Kruskal-Wallis test with two-sided Welch's test as a post hoc analysis and Benjamini-Hochberg FDR correction, which was required due to the non-normal distribution of P values. Q value threshold of 0.05 was assumed. Significance of sample clusters separation was tested with AMOVA on a Bray-Curtis distance matrix derived from the simulated metagenome count table (vegan's vegdist). Functional diversity was calculated as the number of categories at level 0, which is the level of individual Kegg Orthologies (might be understood as functions).
Species diversity was measured as Shannon's diversity (H ′), species richness was measured as observed number of OTUs, and evenness was estimated as Shannon's evennes (J ′). From now on, for the sake of brevity, we will use terms Bdiversity^and Bevenness^in place of their indices names.

Changes of Climatic and Litter Physicochemical Variables over Time
Average daily temperature at the sampling area was below 0°C for the majority of the period preceding t 1 , rose to around 12°C 3 weeks later, and then after a month it fluctuated around 15°C until the end of the experiment. Precipitation was low (1.2 mm daily on average) until day 150; in this period, it was in the form of snow, but the snow cover lasted until day 140. Later, rainfalls were more intensive (3.8 mm daily) until day 230, after which a 3-week period of drought occurred. Final sampling was performed 2 weeks after the onset of drought ( Supplementary Fig. 1).
C/N ratio as well as K level were the highest in t 0 samples and stayed at lower level in older ones ( Table 2). An increasing trend was apparent for Mg, Fe, and Zn, while concentration of Mn was the highest in t 1 and Cu in t 2 samples. pH was similar in all samples; the average was 5.22 ± 0.37. Water content of freshly fallen litter (t 0 samples) was~27%, while directly after removal from the field, t 1 was significantly more humid and contained 53% water and humidity in t 2 was even lower than in t 0 (23%).

Archea Are Absent from Litter Samples
No archaeal reads were recovered during the project. The sequences coming from the libraries prepared using the presumed archaeal primers [37] were classified either as Bacteria (majority as Actinobacteria) or as Eukaryotes (Fungi). As the primers have perfect matches to almost 80% of archaeal sequences in Silva, we think that no Archaea, at least such recoverable with the primer pair used, were present in samples under study.

Addition of Enzymes Does Not Improve DNA Isolation Yield and Microbial Diversity Recovered
Pre-treatment of litter samples with achromopeptidase and lysozyme, chitinase, and all three enzymes together did not affect the copy number of pine 18S rRNA gene in litter metagenomic DNA (P > 0.05, Fig. 1a). The overall DNA yield was not changed by the addition of enzymes (P > 0.05, Fig. 1b). The same applies to the number of copies of bacterial 16S rRNA genes measured with qPCR (P > 0.05, Fig. 1c) and the number of copies of fungal ITS sequences (P > 0.05, Fig. 1d).
Bacterial and eukaryotic diversity (measured as Shannon's H′), species richness as well as evenness were not influenced by modifications of the DNA isolation method, but they differed at various stages of decomposition (Fig. 2). The same applies to community structure, as assessed with AMOVA and ANOSIM performed on Bray-Curtis and Morisita-Horn distance matrices (P > 0.05; data not shown).

DNA Content and Bacterial as well as Fungal SSU Gene Copy Numbers Change During Early Stages of Decomposition
DNA yield significantly increased over time (P < 0.05, Fig. 3a). A rapid decrease of pine 18S rRNA gene copy number in t 1 and t 2 samples was observed (P < 0.05), indicating that most of pine DNA was degraded after 3 months (Fig. 3d). Thus, it seemed that the increase in yield was due to the higher bacterial and fungal cell numbers. Indeed, the copy numbers of both bacteria 16S rRNA genes and fungal ITS sequences, being proxies for bacterial and fungal cell numbers, were greatest in t 2 samples. However, the number of bacterial SSU was significantly higher in t 1 than in t 0 (Fig. 3b), while fungal ITS number was significantly higher in t 2 samples than in t 1 (Fig. 3c).

Microbial Diversity Changes During Litter Decomposition
Bacterial diversity (measured as Shannon's H′) was similar in all samples; it was highest in fresh litter, then dropped after  (Fig. 4a). On the other hand, eukaryotic diversity decreased rapidly over time (Fig. 4a), and the difference between fresh and 8 months old litter was significant. Both observed and estimated total (Chao1) species richness grew with time in the case of bacteria and decreased in the case of eukaryotes (Fig. 4, b, c), while evenness followed the pattern of diversity in both cases (Fig. 4d).

Bacterial and Fungal Community Structure Is Driven by Decomposition Stage, Physicochemical Variables, and Bacteria-Eukaryote Interactions
Regardless of the dissimilarity measure used, bacterial communities coming from individual biological replicates clustered tightly on nMDS plots, and samples at the same stage of decomposition were also located together (Bray-Curtis showed on Fig. 5a). It appeared that the t 1 samples were more  similar to the fresh ones (t 0 ) than to the t 2 samples. On the other hand, eukaryotic communities in t 1 samples were closer to t 2 ones (Bray-Curtis showed on Fig. 5b). ANOSIM as well as AMOVA and PERMANOVA analyses showed that separation of clusters was significant (P < 0.05). PERMDISP test demonstrated that variance was not homogeneous in different sample groups for bacterial community, with highest dispersion found in t 2 samples (P < 0.05), while it was similar in all groups for eukaryotic community (P = 0.848). CCA analysis (Supplementary Fig. 2) identified Zn and Na as significant environmental parameters shaping the bacterial community structure (P < 0.01), while C/N ratio as well as Zn and Mn were factors that significantly influenced eukaryotic community. The measured variables explained 20.1% of inertia (variance) in case of bacteria and 37.2% in case of eukaryotes. Possible influence of eukaryotes on bacterial community and vice versa was assessed with co-correspondence analysis ( Supplementary Fig. 3). Two first COCA axes turned out to be significant in a permutation test (P < 0.05), and they were sufficient (i.e. adding more axes would not significantly increase percent variance explained) according to leave-one-out cross-validation. They explained 22.4% of variance in the bacterial community and 44.15% of variance in the eukaryotic community. Mean distance between bacterial and eukaryotic site coordinates was significantly higher for t 2 samples (0.985, P < 0.05), than for t 0 (0.673) and t 1 (0.646) ones, indicating that stronger influence was exerted in the latter two sample groups.

Bacterial Community in Litter Samples Is Dominated by Proteobacteria
Although bacterial community structure varied greatly in biological replicates, certain general trends were visible. At all levels an increase of rare taxa levels in time was found, which was in line with increasing bacterial diversity. At the level of phylum, the community was dominated by Proteobacteria (68-88%), and an increase of Actinobacteria was visible in t 2 samples (Fig. 6a). Only four phyla (Proteobacteria, Acidobacteria, Bacteroidetes, and Actinobacteria) accounted for the vast majority of reads in all samples. Among classes Alphaproteobacteria dominated in all libraries, and slight decrease of Beta-and Gammaproteobacteria in t 1 and t 2 samples was observed with concomitant increase of Actinobacteria, Sphingobacteriia, and rare taxa levels (Fig. 6b). More changes were visible at the order level, e.g., Rhodospirillales reads were abundant in t 0 reads and significantly less numerous in t 1 and t 2 , Burkholderiales and Xanthomonadales were most frequent in t 1 samples, while Rhizobiales and Sphingobacteriales reads were most frequent in t 2 . The same pattern was observed for families within the abovementioned orders (Fig. 6c) and among genera, where Sphingomonas was the most frequent one, and Pseudomonas, Burkholderia, Granulicella, and Rhizobium were also found in large quantities (Fig. 6d).

Eukaryotic Community in Litter Samples Is Dominated by Fungi
Eukaryotic libraries differed in biological replicates, similarly to bacterial ones. The libraries were dominated by fungal reads (Ascomycota 37-85% and Basidiomycota 7.5-57.8%) with minor quantities of Nematoda (Chromadorea, up to 4.4%) at later stages of decomposition (Fig. 7a). Ascomycotal reads prevailed in the libraries derived from fresh litter, then Basidiomycota were the most numerous in t 1 samples and finally Ascomycota levels increased to~65% in t 2 litter. Each sample type harbored unique eukaryotic community at the genus level (Fig. 7b), in case of t 0 ones hallmarked by Lophodermium and Phoma; Sistotrema, Ceuthospora (Phacidium), Trichoderma, and Athelia being characteristic for t 1 ; and Plectus (nematode), Mycena, as well as Mytilinidion being typical for t 2 . However, due to relatively short 18S rRNA gene fragments being sequenced, only~35 to 80% of reads could be classified down to this level.

Metabolic Capabilities of Bacterial Community Are Different at Different Stages of Decomposition
PICRUSt analysis was performed to reveal potential functions encoded in genomes of bacteria whose 16S rRNA gene fragments were sequenced in our study. Nearest Sequenced Taxon Index (NSTI) analysis indicated that in most cases, a sequenced genome from the same genus could be found (NSTI < 0.05).
Overall diversity of functions was lower in the t 2 samples than in t 1 and t 0 (P < 0.01; Fig. 8a). PCA analysis showed that functional composition of metagenomes was different in t 0 , t 1 , and t 2 ( Fig. 8b; P < 0.05, AMOVA); however, this effect was partially due to differences in variance (P < 0.01, PERMDISP). The pathways whose shares differed significantly between various sample types belonged mainly to BMetabolism^supercategory, but there were also BSignaling Molecules and Interaction^as well as BTransport and Catabolism^ (Fig. 9).
Carbohydrate, amino acid, energy, lipid, xenobiotic, and terpenoid metabolism were the most frequent categories differing significantly. Genes encoding proteins participating in carbohydrates, terpenoids and xenobiotic metabolism were predicted to be more frequent in t 2 communities than in the other two, while those engaged in lipid metabolism were estimated to be more frequent in t 2 than in t 1 only. Only the genes involved in energy metabolism were more probable to be encoded by organisms forming t 0 community. A B 66 Gołębiewski M. et al. At the level of individual genes, these involved in various sugar transport (PTS system) and certain engaged in metabolism (e.g., sugar kinases, epimerases) as well as regulation (sugar utilization operon regulatory proteins) were predicted to be most abundantly represented in genomes of organisms from t 0 samples, virtually absent from t 2 samples with t 1 in the middle. The same applies to other genes involved in import and utilization of soluble substrates such as amino acids, lipids, or amines, while β-glucosidases, potentially engaged in cellulose degradation, were predicted to be most frequent in t 1 samples. In contrast, genes involved in biosynthesis of vitamins and cofactors, degradation of aromatic compounds (oxygenases and oxidases responsible for degradation of sterols and phenolic compounds such as catechol, 2,4dichlorophenoxyacetate or 4-hydroxyacetophenone), as well as Fe 3+ transporting ATPases were predicted to be most frequent in genomes of organisms dwelling in t 2 samples.

Discussion
The short time span of our experiment together with the small number of sampling time points causes our results to be rather a short sequence of snapshots than a time course of microbial succession on litter. Therefore, we avoid interpreting the results as Btrends^and concentrate on differences between pairs of time points. Nevertheless, we think that the results provide an interesting insight into early stages of pine litter decomposition at the molecular level. We observed no influence of enzymatic treatment neither on the amount of DNA isolated nor on microbial diversity. This might have been caused by the effectiveness of bead beating, making the action of enzymes pointless. In our previous study, we observed the influence of enzymes, as the method of isolation was based on less effective thermal/ chemical lysis [61]. However, an alternative explanation is possible, assuming that cell envelopes of microbes in our samples were not digested by the enzymes used. Albeit possible, it is unlikely at least in the case of lysozyme, taking into consideration high numbers of Proteobacterial reads generated from the samples. We found that the enzyme worked well on E. coli DH10B in the buffer used for digestion. However, no increase in fungal DNA yield could be an effect of digestion conditions that were suboptimal for the chitinase used (according to the producer the optimal pH is 4, while the buffer used for digestion had pH 8). It seems that there is no reason for enzymatic pre-treatment of litter samples prior to beadbeating-based DNA isolation. Other studies report significant influence of lysozyme digestion on soil DNA yield, even in combination with bead-beating [62]; however, the beadbeating method was not as effective as in contemporary kits using dedicated bead-beaters.
DNA content of litter samples significantly increased during the studied decomposition period. This increase was caused by greater numbers of bacterial and fungal cells, as we found that pine DNA was almost completely degraded in t 1 samples. This was confirmed by results of qPCR analysis of  bacterial 16S rRNA and fungal ITS copy numbers. The number of bacterial sequences increased first (in the t 1 samples), supporting our hypothesis that fresh litter is first colonized by bacteria, albeit, alternatively, this increase might have been caused by growth of phyllospheric organisms. It seems plausible that bacteria are responsible for pine DNA degradation. We have found no archaeal reads in our libraries, which might mean that there were no Archaea in our samples. This is supported by results of a metaproteomic study concerning beech litter decomposition [16]. Nevertheless, it is possible that Archaea in pine litter might belong to a group whose 16S rRNA sequences differ from those amplifiable with the primer system we used, e.g., Nanoarchaeota, Altiarchaeota, or Diapherotrites or, as in all phyla the coverage was well below 90%, they could be members of generally amplifiable groups having non-amplifiable variants. There are reports of Archaea being cultured from Scots pine ectomycorrhizas (Methanolobus, Halobacterium, and unknown member of 11.c Crenarchaeota) [63], and~2% of sequenced transcripts from spruce litter were coming from Archaea [26] and both facts support the latter possibility.
Shannon's diversity, evenness, and species richness were consistently higher for bacterial communities than for fungal ones, which might be caused by the number of bacterial taxa being around one order of magnitude greater than the number of fungal ones in soils [64]. Moreover, bacterial diversity grew over time, while fungal communities became less diverse. This effect might have been caused by three factors: litter humidity, temperature, and time (decomposition stage).
Although precipitation was lower in the t 0 -t 1 period than in t 1 -t 2 one, water content was the highest in t 1 samples, which was probably caused by low temperature. During the initial period of the experiment (t 0 -t 1 ), temperature was~0°C with occasional snow, and diversity did not change significantly, neither in the case of fungi, nor bacteria, which suggests weak selection. Later, warm period was hallmarked by bacterial species richness increase, showing possible colonization of litter by bacteria from the surrounding environment. Evenness was lower in t 2 samples causing diversity to stay at the t 0 level. In the case of eukaryotic community diversity, species richness and evenness decreased significantly in the t 2 samples, suggesting strong selection.
t 0 samples harbored a bacterial community similar to pine phyllosphere community [65], with high shares of Proteobacterial sequences, particularly members of the Acetobacteraceae and Pseudomonadaceae families being the hallmarks of phyllosphere. Differences might have resulted from tree species being different (Scots pine vs. limber pine), geographical distance (Central Europe vs. California), and physiological state of needles (dry, fallen vs. fresh ones). The eukaryotic t 0 community consisted of Fungi, out of which Ascomycotal class Pezizomycotina prevailed, which was similar to fresh pine needles community studied by Millberg and colleagues [66] and unlike in senescent oak leaves, where Dothideomycetes prevailed [17]. Again, the differences might have resulted from geographical distance (Sweden vs. Poland) and physiological state of the needles. Eukaryotic sequences that could be identified down to the genus level belonged mainly to known phyllospheric fungi, such as Phoma, Trichoderma, or Lophodermium, and the latter might act as an early decomposer in litter [67,68]. Low number of sequences was assigned to Naemacyclus (Helotiales, synonym Cyclaneusma), a genus comprising fungi causing needle cast [69]; thus, it is plausible that certain amount of collected needles were shed prematurely due to its action. Differences between t 1 and t 0 samples were profound, in spite of temperature being below 0°C most of the time. After 3 months of field incubation, the samples harbored communities that were drastically different from t 0 ones; this fact demonstrates that colonization of fresh litter by soil and older litter-inhabiting organisms is rapid, regardless of harsh environmental conditions. This is in line with results of earlier studies on mass loss, e.g., [70] and with reports concerning decomposition of broadleaf litter [17,19,20]. Reads coming from phyllospheric organisms were less abundant than in the t 0 samples (e.g., members of Acetobacteraceae, as well as fungi Lophodermium and Phoma), showing that, in spite of the lack of diversity decrease, selection operates at this stage of litter decomposition, but organisms unable to survive are replaced by colonizers. Obvious fungal colonizers were members of Sistotrema, ectomycorrhizal/saprophytic fungi of the Cantarellales order (Basidiomycota, Agaricomycetes), and Ceuthospora (synonym Phacidium, of Helotiales (Ascomycota, Dothideomycetes), [71], a genus comprising plant pathogenic fungi causing various diseases from needles rot in conifers to fruit rot in cranberries, apples, or pears [72,73]. Reads derived from these organisms were virtually absent from t 0 samples and constituted over 50% of all eukaryotic reads in t 1 ones. Small numbers of reads obtained from t 0 samples might be explained by spores sedimenting from air on the collected needles. The eukaryotic community at this stage of decomposition was more similar to the initial one than the bacterial one, suggesting that bacteria were able to grow more actively under winter conditions. C:N ratio was highest in t 0 samples, indicating nitrogen depletion, probably due to organic nitrogen being consumed by microbes.
To no surprise, at the end of the experiment, the samples harbored communities different from both t 0 and t 1 ones. As seasonal differences in microbial community composition in spruce litter were not found to be large [26], we suggest that the differences visible in our data should be attributed to decomposition stage. Phyllospheric organisms were even less abundant in t 2 than in t 1 , typical soil organisms could be found, such as nematodes of the Chromadorea class or bacteria of the Rhizobium genus, and the litter was colonized by organisms belonging to Dothideomycetes (Fungi, Ascomycota) and Sphingomonadaceae as well Bradyrhizobiaceae (Bacteria, Alphaproteobacteria), probably capable of degradation of recalcitrant compounds (lignocellulose, lignin) prevailing at this stage of decomposition [74]. The presence of nematodes might be one of the causes of slower increase of bacterial 16S rDNA counts, as the most abundant phylotype of this group (Plectus) is a bacterivore. This is supported by the results of cocorrespondence analysis, as the influence of eukaryotic community on the bacterial one is significant and explains~10% of variation. On the other hand, interactions among Eukaryota cannot be excluded, e.g., certain fungivorous nematodes most abundant in the t 1 samples (e.g., members of the Aphelenchoides genus [75]) might be involved in much lower Sistotrema abundance in the t 2 samples. The influence of bacteria on the eukaryotic community was even stronger (~22%). It might be explained, e.g., by providing nitrogen (by N 2 fixation) that seems to be a limiting factor, particularly in t 0 samples. This is confirmed by the results of CCA, wherein C/N was identified as a significant environmental variable. It is also possible that bacteria provided fungi with phosphorus, another limiting nutrient, by solubilizing its soil resources [76]. As the measured environmental variables were responsible for explaining of over 37% of variance, it seems plausible that they were the key drivers of microbial community structure changes. Bacterial community structure at the level of genus was surprisingly similar to this found at the corresponding stages of decomposition in oak leaf litter [77], where Pseudomonas, Sphingomonas, and Burkholderia were among the dominating genera. Differences were mainly quantitative and were more pronounced in the phyllospheric communities. The most remarkable difference was the lack of Duganella (Oxalobacteriaceae, Betaproteobacteria) and Frigoribacterium (Micrococcaceae, Actinobacteria) in pine needles; these genera were replaced by unclassified members of the Acetobacteraceae family. Subsequent stages of pine litter decomposition harbored less Pedobacter. This fact suggests that, while plants assemble unique phyllospheric communities, decomposers' assemblages, at least bacterial, may be more generic and similar regardless of the quality of litter.
Fungal communities decomposing pine litter are less similar to those degrading oak leaf litter described by Voriskova and colleagues [17]; however, general picture looks similar: phyllospheric communities are rapidly replaced by distinct communities characteristic for particular stages of decomposition. Certain fungi are abundant both in oak and pine decomposing litter, e.g., Athelia, Rhodotorula, and Sistotrema. Interestingly, Sistotrema, most abundant in 12-month-old oak litter [17], displays abundance peak in t 1 samples. Thus, it is possible that fungal succession on decomposing pine litter Bovertakest he one on oak leaves.
A model of litter decomposition assumes that soluble compounds are degraded first, then degradation of hemicelluloses follows and the last stages comprise degradation of celluloselignin complex leading to increasing lignin/cellulose ratio [78]. We used PICRUSt to model possible metagenomes of decomposing litter to learn if the probable gene content supports this model. However, one should bear in mind that the results obtained with PICRUSt are only an approximation of the Btruem etagenome, and many factors, such as all biases influencing the underlying 16S rRNA gene fragment sequencing or differences in gene content of closely related organisms due to lateral gene transfer, may skew the results. Moreover, the set of functions derived from PICRUSt is only potential one, i.e., we do not know if the functions are actually expressed. Another limitation is that, due to the lack of a eukaryotic database, it was not possible to perform such an analysis for eukaryotic sequences. Therefore, the results should be treated with caution.
We expected a shift in community composition from organisms capable of various carbohydrates utilization to cellulose-degrading specialists, to lignin degraders, although the latter only to some extent, as our experiment concentrated on the first 8 months of litter decomposition. Our results support this view, as cellulolytic genes (β-glucosidases) were most abundant in t 2 samples, and certain genes that might be involved in lignin degradation (oxygenases and oxidases) were also found. The genes engaged in active transport of cellobiose (cellulose degradation product) into the cell (PTS system) were most frequent in t 0 samples and permeases allowing for passive import of this sugar were most numerous in t 1 samples. This picture might mean that during the initial stage of litter decomposition (t 0 ) after depletion of other easily accessible sugars, active transfer of cellobiose is required, most probably due to the low concentration of this compound; then, fungi start producing cellulases causing increase of the cellobiose concentration, which allows passive transfer; and finally, bacterial cellulolytic activity may supplement fungal activity (t 2 ). The latter claim was partially supported by results of Zifcakova and colleagues, who found that Fungi were responsible for only half of carbohydrates metabolism-involved transcripts production in spruce litter [26]. Given the importance of coniferous litter decomposition and taking into account all the limitations of the methodology as well as the general scarcity of molecular studies devoted to coniferous litter decomposition, more work is needed to understand which organisms are responsible for particular processes during decomposition.

Conclusions
There is no need to supplement contemporary bead beatingbased kits with cell wall-degrading enzymes. Bacterial and fungal diversity changes differently during initial stages of decomposition. Microbial succession on decomposing pine litter is rapid and initial phyllospheric communities are replaced with distinct assemblage characteristic for particular stages of decomposition. Changes in communities seem to be driven not only by physicochemical variables of litter but also by interactions among bacteria, fungi, and nematodes.
Author Contributions Conceived the study-MG, MN, designed the experiments-MG, contributed financial support and reagents-MN, MG, ATa, ATr, executed the experiments-MS, ATa, ED-S, MG, performed bioinformatic and statistical analyses-MG, interpreted data-MG, ATa, MN, wrote the draft of the manuscript-MG. All of the authors revised the manuscript and approved its content.
Funding Information This study was performed within an OPUS grant of the Polish National Science Center (No. UMO-2012/05/NZ8/001362) and supported by DS758, DS759, and DSC funds of the Institute of Environmental Sciences, Jagiellonian University and by statutory funding from Nicolaus Copernicus University (to MG and ATr).
Open Access This article is distributed under the terms of the Creative Comm ons Attribution 4.0 International License (http:// creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. abundance. Dot and whiskers plots to the right of the barplots show 95% confidence intervals on the difference between respective samples, dot location indicates the mean, and the color of a given dot shows in which sample the proportion was higher. Bonferroni-corrected P value is given to the right of dot and whiskers plots