Recent years with growing availability of sequencing technologies and respective instrumentation have seen increasingly expanding investigation of bacterial community diversity in various environments. Yet little is known about the structure and richness of bacterial communities in forest soil, especially how the main forest-forming tree species can shape microbial community composition and functioning, and to what extent soil microorganisms found in the root zone can influence plant growth and development. Rhizosphere, which is a small volume of soil influenced by root exudates, is the main interface between plant and soil. Root exudates are the main source of carbon substrates available for soil microorganisms [Kruse et al., 2013], supporting dozens of billions bacterial cells per gram of root [Egamberdieva et al., 2008] and tens of thousands of prokaryotic species [Mendes et al., 2011]. Even small shifts in root exudates composition may bring about drastic changes in microbial community composition [Falcini et al., 2003; Badri et al., 2009], which is currently even thought to be the main mechanism deployed by plants to control pathogen development in their rhizosphere [Berendsen et al., 2012]. Tree species were shown to impose a selective force on the composition of the rhizosphere bacterial communities [Gomes et al., 2014], the phenomenon appearing to be plant species specific. Thus it is logical to expect that differences in soil microbial community composition and structure may be more pronounced in the rhizosphere. So our study aimed to assess the difference in rhizosphere bacterial community (16S rRNA-based) structural diversity between two pine species (Pinus sibirica Du Tour and Pinus koraiensis Sieboldd et Zucc.), growing in the long-term field provenance experiment in the south of Krasnoyarsk region (Russia). Materials and methods Experimental setup and soil sampling The provenance experiment with different pine species and their climatypes was set up in the Ermakovo Forestry (53°10´ N, 95° 20´ E) of the Krasnoyarsk region, Russia, in 1983. The experiment is described in details elsewhere [Naumova et al., 2014]. Briefly, the Siberian pine (Pinus sibirica, PS) was represented by 3 climatypes, namely, from the Kemerovo region (K), Tomsk region (T) and Ermakovo (E) of the Krasnoyarsk region. The latter climatype commonly grows in the immediate vicinity of the experiment plot and by default is better adapted to the environment. The bulk samples collected from under these three climatypes statistically represented three individual replicates to characterize Siberian pine rhizosphere. The Korean pine (Pinus koraiensis, PK) was represented by 2 climatypes, namely, from the Khabarovsk region (Kh) and from the Primorsk region (P) and by default has been growing in the non-native (much drier and colder) environment. The bulk samples collected from under these two climatypes statistically represented two individual replicates to characterize Korean pine rhizosphere. Pine trees’ characteristics are shown in Table 1. Table 1 Some tree growth and development properties of the two pine species (Siberian pine and Korean pine), mean ± s. e. m. Property Siberian pine Korean pine Tree height, m 8.30 ± 0.65 9.85 ± 1.35 Trunk diameter, cm 11.7 ± 0.2 12.8 ± 0.2 Tree height yearly increment, cm 43 ± 2 48 ± 0 Crone diameter, m 2.5 ± 0.1 3.4 ± 0.4 Number of branches per verticil 7.3 ± 0.4 4.5 ± 0.3 Survival 53 ± 16 81 ± 2 Rhizosphere soil was sampled in August 2013 by collecting roots from the six trees per plot (0-20 cm layer, i.e. immediately below the litter layer and with maximal root density) at the distance of 60 cm from the tree row within the crone zone. Rhizosphere soil was collected as the soil remaining on 1-3 mm pine roots after their gentle shaking [Zhao et al., 2010]. One composite sample was mixed from roots of the 6 individual soil monoliths from each plot. After sampling soil was brought into the laboratory for DNA extraction and chemical analyses. Some soil chemical and microbiological properties, described in detail earlier [Makarikova et al., 2014; Naumova et al., 2014], are shown in Table 2. Table 2 Some chemical and microbiological properties of the rhizosphere soil of the two pine species (Siberian pine and Korean pine), mean ± s. e. m. Property Siberian pine Korean pine pH 5.96 ± 0.05 5.65 ± 0.64 SOC* 5.4 ± 0.4 4.6 ± 0.8 DOC, µg С g-1 soil 97 ± 17 80 ± 23 SUVA, l∙ cm-1∙ g-1 DOC 132 ± 20 99 ± 9 STN, 0.140 ± 0.001 0.099 ± 0.005 C/N 33 ± 3 45 ± 7 NO3-, mg N·kg-1 soil 2.6 ± 0.7 1.5 ± 0.2 NO2-, mg N·kg-1 soil 0.38 ± 0.03 0.21 ± 0.03 NH4+, mg N·kg-1 soil 4.3 ± 0.9 2.2 ± 0.8 P2О5, mg·kg-1 soil 5.1 ± 0.6 3.8 ± 0.8 SMBC, μg С · g-1 o.d. soil 590 ± 89 413 ± 76 SMBC/SOC, 1.1± 0.1 0.9 ± 0.0 SIR, μl · hr-1 · g-1 o.d. soil 13.9 ± 1.4 10.3 ± 0.9 QR 0.33 ± 0.03 0.40 ± 0.13 Qmet, μg С-СО2 ·mg SMBC-1·hr-1 3.8 ± 0.1 4.8 ± 1.2 *Abbreviations used: SOC - soil organic C, DOC - dissolved organic carbon, SUVA - specific UV absorbance, STN - soil total nitrogen, C/N - soil organic matter C/N ratio, SMBC - soil microbial biomass C, SMBC/SOC - proportion of soil microbial biomass C in the soil organic C, SIR - soil substrate-induced respiration, QR - respiratory coefficient, Qmet - metabolic coefficient. Extraction of total nucleic acid from soil Total DNA from 0.40 g of soil was extracted using the DNA isolation Kit (DNA Spin Kit for Soil™, MO Bio laboratories, Inc., USA) as per manufacturer’s instructions. The bead-beating was performed using mini bead-beater (Stanford, USA), for 45 s at 5000 rpm. No further purification of the DNA was needed. The quality of the DNA was assessed using agarose gel electrophoresis. 16S metagenomic sequencing The V3-V4 region of the 16S rRNA genes was amplified with the primer pair 343F (5’-CTCCTACGGRRSGCAGCAG-3’) and 806R (5’-GGACTACNVGGGTWTCTAAT-3’) combined with Illumina adapter sequences, a pad and a linker of two bases, as well as barcodes on the primers [Caporaso et al., 2012]. PCR amplification was performed in 50 μl reactions containing 0.7 U Phusion Hot Start II High-Fidelity and 1× Phusion GC buffer (Thermo Fisher Scientific), 0.2 μM of each forward and reverse primers, 10 ng template DNA, 2.3 mM MgCl2 (Sigma-Aldrich) and 0.2 mM of each dNTP (Life Technologies). Thermal cycling conditions were as follows: initial denaturation at 98°C for 1 min, followed by 30 cycles of 98 °C for 15 s, 62 °C for 15 s, and 72 °C for 15 s, with final extension at 72 °C for 10 min. A total of 200 ng PCR product from each sample was pooled together and purified through MinElute Gel Extraction Kit (Qiagen). Sample libraries for sequencing were prepared according to the MiSeq Reagent Kit Preparation Guide (Illumina) and the protocol described previously [Caporaso et al., 2011; 2012]. Sample denaturation was performed by mixing 4.5 µl of combined PCR products (4 nM) and 4.5 µl 0.2 M NaOH. Denatured DNA was diluted to 14 pM and 510 μl mixed with 90 μl of 14 pM Phix library. A total of 600 µl sample mixture, together with customized sequencing primers for forward, reverse, and index reads, were loaded into the corresponding wells on the reagent cartridge of the 500-cycle PE kit and run for 2x250 bp paired-ends sequencing on MiSeq Illumina sequencer at the SB RAS Genomics Core Facility (ICBFM SB RAS, Novosibirsk, Russia). Bioinformatic analysis 300 PEs were overlapped, quality filtered (QV>20) and demultiplexed. Any overlapped reads with ambiguous sites and also without primers were removed. All steps above were done using CLC GW 7.0 (CLC Bio). RDP Classifier 2.10 was used for taxonomic analysis [Wang et al., 2007]. Before analysis chimera checking was performed by usearch 7.0 [Edgar et al., 2011]. Statistical analyses (PCA, ANOVA) were performed by using Statistica v.6.1 software package. Bacterial genera sequences diversity in the rhizosphere of the studied pines was compared by calculating diversity indices using PAST software v.2.17. Results and Discussion 16S rRNA-based bacterial diversity Analysis of 16S rRNA-based composition and structural diversity of bacterial communities found 196 ± 62 genera in PS rhizosphere and 168 ± 3 genera in PK rhizosphere. The number of genera contributing at least 1% into the total number of sequences was 20 under PS and 17 under PK. So the overwhelming majority, i.e. about 90 of the genera sequences, were minor members of the rhizosphere communities, contributing mostly less than 1% into total community sequences. The dominant genera sequences were similar under both pine species. The list of thus defined dominants, overall 22 genera with their sequence percentage relative to the total sequences exceeding 1%, is shown in Tab.3. Both rhizosphere soils were dominated by Bradyrhizobium ssp. The Bradyrhizobium genus, a member of Alphaproteobacteria, has been well known as nitrogen-fixing, leguminous roots nodulating bacteria of agricultural importance [Rodríguez-Echeverría et al., 2014]. We found few researches reporting the natural occurrence of this bacterial genus in soil under trees: for example, in natural and artificial forest stands in Brazil [Pereira et al., 2008]. Apparently, in our experiment carbon into the pine rhizosphere was enough to sustain Bradyrhizobium ssp. to perform nitrogen fixation. The number of Bradyrhizobium sequence reads in PK rhizosphere soil was 1.3 times higher as compared to the PS one; if the ratio of SSU 16S rRNA gene copies to nif gene copies and activity was similar in both rhizosphere communities studied, it might suggest higher nitrogen fixation rate in PK rhizosphere. However, soil nitrogen contents, both mineral and organic, were found to be higher in PS rhizosphere soil (Table 2). The discrepancy may be explained by the higher presence of other diazotrophs among the PS rhizosphere sequences, as well as by higher nitrogen uptake and phytomass sequestration by more actively growing PK plants. The Spartobacteria (genus of uncertain taxonomic placement, second dominant under PS and forth dominant under PK) belongs to Veruccomicrobia phylum that is widely distributed in aquatic ecosystems, soils and rhizospheres [Buckley & Schmidt, 2001; Sangwan et al., 2005]. Despite increasing evidence that they are a numerically abundant component of soil microbial communities, very little is known about the ecological significance of Verrucomicrobia, as its culturable members were isolated from plant-soil environment not long ago [Nunes da Rocha et al., 2009], and most of the ecological information about this phylum is based on Verrucomicrobia 16S rRNA gene sequence distribution over different ecosystems. Recently it was shown that Veruccomicrobia phylotypes, dominating bacterial community in sediments of an arctic fjord, were responsible for hydrolysis of polysaccharide substrates at high rates [Cardman et al., 2014]. So the dominance of the genus in the pine rhizosphere soils, found in this study, may be related to the rhizodeposition of oligo- and/or polysaccharides. Earlier in field experiments the presence of plants was shown to favour the development of a range of bacterial phyla, including Veruccomicrobia [Nunes da Rocha et al., 2013]. Table 3 The dominant bacterial genera in the rhizosphere soil of the two pine species (Siberian pine and Korean pine) as judged by the percentage of the specific sequences in the total number of sequences (mean ± s. e. m.) Genus Siberian pine Korean pine Bradyrhizobium 10.5 ± 0.8 14.1 ± 1.7 Acidobacteria (subdivision 2)* 7.9 ± 0.8 13.0 ± 0.6 Acidobacteria (subdivision 1) 6.2 ± 1.6 9.1 ± 2.0 Spartobacteria (gis)** 10.1 ± 1.4 7.5 ± 0.6 Gaiella 8.1 ± 2.8 6.7 ± 0.6 Acidobacteria (subdivision 3) 3.5 ± 1.3 6.0 ± 1.7 Mycobacterium 2.6 ± 0.5 4.5 ± 0.6 Phenylobacterium 3.2 ± 0.3 3.3 ± 0.1 Mucilaginibacter 2.9 ± 0.7 3.2 ± 1.6 Acidobacteria (subdivision 6) 8.1 ± 0.7 3.1 ± 1.0 Conexibacter 1.9 ± 0.2 2.8 ± 0.2 Burkholderia 0.8 ± 0.6 2.6 ± 0.3 Ktedonobacter 1.6 ± 0.3 2.5 ± 0.2 Candidatus Solibacter 1.1 ± 0.3 1.9 ± 0.2 Granulicella 0.7 ± 0.3 1.8 ± 0.3 Subdivision 3*** (gis) 1.7 ± 0.3 1.7 ± 1.3 Gemmatimonas 1.4 ± 0.2 1.1 ± 0.0 Acidobacteria (subdivision 16) 2.5 ± 0.1 0.7 ± 0.2 Sphingomonas 1.2 ± 0.4 0.5 ± 0.2 Acidobacteria (subdivision 4) 1.3 ± 0.1 0.5 ± 0.2 Phaselicystis 1.0 ± 0.1 0.4 ± 0.0 Ilumatobacter 1.5 ± 0.0 0.1 ± 0.1 Total 79.9 87.2 * The subdivision number is shown in brackets; ** genus incertae sedis; *** - of Veruccomicrobia. The subdivision 2 Acidobacteria genus is second dominant genus in PK rhizosphere (Tab.3). The Acidobacteria phylum is highly diverse, physiologically active and well represented in soil environment [Naether et al., 2012]. Like Veruccomicrobia, members of this phylum have been difficult to culture, so the knowledge of their metabolic capabilities and ecological potential is yet very limited. However, some data from whole-genome analysis suggested that at least some Acidobacteria are versatile heterotrophs capable of using both plant polymers and readily oxidizable carbon [Eichorst et al., 2011; Rawat et al., 2012]. Some of its species were found to harbor increased numbers of mobile element genes compared with the other Acidobacteria genomes [Challacombe & Kuske, 2012]. This could be due to a particular need for enhanced competitive ability under extremes of moisture, temperature, geochemical conditions and other environmental conditions. ANOVA analysis of genera sequence abundance revealed that out of the major bacterial community components pine species affected Acidobacteria subdivision 16 (P ≤ 0.045), as well as several minor components (Acidobacteria sundivision 10, Actinomadura, Marmoricola, Iamia and Methylocella). Their ecological significance should be studied more closely as, albeit minor, they still can perform key ecological functions. Based on genera sequences, the PS rhizosphere bacterial community appeared to be richer and more even as compared to the PK rhizosphere soil (Tab.4): the total numbers of genera, as well as genera richness indices (Shannon, Menhinick, Margalef, Fisher’s alpha) were higher in PS rhizosphere, i.e. rhizosphere of pine species growing in its native environment. As in many situations the Fisher’s alpha index is approximately equal to the number of taxa represented by a single individual, the minor fraction of the PS rhizosphere bacterial sequences was more diverse as compared to the PK one. The dominance indices (dominance, Berger-Parker) were higher, while evenness indices (evenness, Simpson, equitability), were lower in the PK rhizosphere soil. Thus bacterial DNA diversity comparison allows concluding that plant species introduced into a novel habitat can shape soil bacterial community structure differently as compared to the indigenous species. Table 4 Bacterial community diversity indices of the rhizosphere soil of the two pine species (Siberian pine and Korean pine) Index Siberian pine Korean pine p (bootstrap) p (permutation) Number of genera 278 206 0.001 0.001 Number of sequences 13714 15578 0 0 Shannon 3.52 3.21 0.001 0.001 Dominance 0.05 0.07 0.001 0.001 Simpson 0.95 0.93 0.001 0.001 Evenness 0.12 0.12 0.737 0.753 Menhinick 2.4 1.7 0.001 0.001 Margalef 29.1 21.2 0.001 0.001 Equitability 0.63 0.60 0.001 0.001 Fisher’s alpha 49 34 0.001 0.001 Berger-Parker 0.11 0.15 0.001 0.001 The principal components analysis performed on the matrix with diversity indices as active variables, soil microbial biomass and activity characteristics as supplementary variables, i.e. variables not used explicitly for analysis, but for interpretation purposes, and rhizosphere samples as objects (Fig.1) grouped pine species separately and showed that the dominance indices closely correlate with glucose-induced respiration and metabolic quotient of soil microbial biomass, both known to be higher stressed conditions, while the evenness indices correlate with the soil microbial biomass carbon. Interestingly, the total number of sequences correlated positively and very closely with soil respiratory activity. Based on covariance [Abdi & Williams, 2010] principal components analysis of the data matrix with bacterial genera sequence numbers as variables and rhizosphere soil samples as objects allows to conclude that increased abundance of such bacterial genera sequences as Acidobacteria (subdivisions 16 and 6), Gaiella and a Spartobacteria genus of uncertain placement were mostly associated with PS rhizosphere, while Acidobacteria (1, 2 and 3), Bradyrhizobium and Mycobacterium were more closely related to PK rhizosphere (Figure 2). Overall, in the studied rhizosphere soil samples more than 90 of the sequences were found to belong to 11 classes, with more than 50 of the sequences representing just three classes, namely Alphaproteobacteria, Actinobacteria and Spartobacteria. As for the phyla, 46% of the total sequences belonged to Proteobacteria, while 21% and 17% represented Acidobacteria and Actinobacteria, respectively. Proteobacteria phylum has been well known to be the key component of soil bacterial communities under coniferous forests [Sun et al., 2014], and high profile of Acidobacteria in forest bulk and rhizosphere soils has been increasingly revealed recently [Mukherjee et al., 2013; Shi et al., 2015]. Figure 1. Projection of soil microbiological characteristics and diversity indices onto the plane of the first two principal components. Active variables: Gen - total number of genera, Seq - total number of sequences, Mar - Margalef index, Fa - Fisher’s index, Sha - Shannon index, Sim - Simpson index, Equ - equitability, Men - Menhinick index, Eve - evenness, Dom - dominance, BP - Berger-Parker index. Supplementary variables: SMBC - soil microbial biomass C, SMBC/SOC - proportion of soil microbial biomass C in the soil organic C, SMBN - soil microbial biomass N, SMBN/STN - proportion of soil microbial biomass N in soil total N, SMB(C/N) - soil microbial biomass C/N ratio, CO2 - soil basal respiration, SIR - soil substrate-induced respiration, QR - respiratory coefficient, Qmet - metabolic coefficient. Siberian pine climatypes: E, Ermakovo; K, Kemerovo; T, Tomsk; Korean pine climatypes: Kh, Khabarovsk, P, Primorsk. Figure 2. Principal components analysis of the number of SSU 16S rRNA gene sequences of different bacterial genera: location of bacterial genera and rhizosphere soil samples in the plane of the first two principal components (PC 1 and PC2, factor loadings). Siberian pine climatypes: E, Ermakovo; K, Kemerovo; T, Tomsk; Korean pine climatypes: Kh, Khabarovsk, P, Primorsk. Table 5 Factor (pine species) contribution into the total variance of SSU 16S rRNA gene sequence abundance in the rhizosphere soil of the two pine species (Siberian pine and Korean pine) and its significance (p) as analysed by ANOVA Class Factor (pine species) Significance (p) Class Factor (pine species) Significance (p) Acidobacteria (10)* 89 0.02 Cytophagia 29 0.35 Acidobacteria (16) 79 0.05 Ktedonobacteria 22 0.42 Anaerolineae 66 0.10 Caldilineae 22 0.42 Acidobacteria (3) 63 0.11 Planctomycetia 22 0.42 Negativicutes 63 0.11 Acidobacteria (17) 20 0.46 Acidobacteria (5) 61 0.12 Verrucomicrobiae 17 0.48 Acidobacteria (4) 60 0.12 Spirochaetia 17 0.50 Thermoleophilia 58 0.13 Acidobacteria (18) 17 0.50 Flavobacteriia 56 0.14 Holophagae 17 0.50 Thermomicrobia 56 0.15 Fimbriimonadia 17 0.50 Deltaproteobacteria 51 0.18 Acidobacteria (11) 13 0.55 Acidobacteria (2) 50 0.18 Bacteroidetes ** 12 0.57 Acidobacteria (13) 50 0.18 Subdivision3 11 0.59 Acidobacteria (7) 50 0.18 Cyanobacteria 11 0.59 Acidobacteria (1) 48 0.20 Betaproteobacteria 9 0.63 Acidobacteria (6) 42 0.24 Spartobacteria 9 0.63 Chlamydiia 41 0.24 Acidobacteria (15) 7 0.67 Fusobacteriia 38 0.27 Gammaproteobacteria 6 0.70 Acidobacteria (12) 38 0.27 Opitutae 3 0.78 Acidobacteria (22) 37 0.28 Sphingobacteriia 3 0.80 Nitrospira 34 0.31 Chthonomonadetes 2 0.81 Acidobacteria (25) 31 0.33 Gemmatimonadetes 1 0.85 Bacilli 31 0.33 Actinobacteria 0 0.93 Alphaproteobacteria 29 0,35 Clostridia 0 1,00 * The subdivision number is shown in brackets; ** genus incertae sedis Bacterial community and pine growth and development The variance of SSU 16S rRNA gene sequence abundance in the rhizosphere soil associated with pine interspecies heterogeneity accounted for ≥ 50% for 14 classes (Tab.5), of which 8 classes were those of Acidobacteria. The difference between pine species was found to be statistically significant for Acidobacteria subdivisions 10 and 16. Soil bacteria dwelling in the rhizosphere impact plant growth and development. As coniferous tree species are increasingly used for reforesting, remediation and designing purposes, it will be helpful to understand the benefits of bacteria for tree fitness. To obtain an idea of such relationship, we performed principal components analysis of the data matrix with tree growth and development parameters (tree height and its annual increment, trunk and crone diameters, and the number of branches per a vertical, as in Tab.1) as variables for analysis and principal components used to obtain Fig.1, i.e. based on SSU 16S rRNA gene sequence abundance, as supplementary variables, visualized the structure of relationships among the variables, on the one hand, and tree climatypes and species, on the other hand (Fig.3). Judging by tree height, its annual increment, trunk and crone diameters, Korean pine, despite the non-native soil and climatic environment, has been growing much better as compared to Siberian pine (Table 1). Pine height and annual height increment is positively correlated with sequence-based PC2, i.e. with such bacterial genera as associated with the positive pole of the PC2 in Fig.2, namely Blastococcus, Enhydrobacter, Mycobacterium and Acidothermus. Tree trunk and crone diameter was negatively correlated with sequence-based PC1, i.e. with such bacterial genera as associated with the negative of the PC1 in Fig.2, namely Acidobacteria 1, 2, 3 and Bradyrhizobium. However, we can not draw any conclusions about the cause-and-effect nature of this relationship as both a) difference in phytomass chemistry and amount can influence soil bacterial community composition and structure [Falcini et al., 2003; Berendsen et al., 2012], b) the latter, reciprocally, by changing soil environment, e.g. organic matter mineralization and plant nutrients availability, as well as abundance and/or activity of other microbial species, can affect tree growth and development, and c) soil fungi may be of more importance in this aspect. In this study we could not draw any conclusions about the relationship between pine intraspecies variation and bacterial DNA diversity in their rhizosphere since soils from under different climatypes were not individually replicated; however, one can say that bacterial diversity variability associated with climatypes seemed comparable with that associated with pine species. It should be noted that the results presented here, are based on analysis of 16S rRNA gene sequence abundance in the rhizosphere soil of two pine species and can not be directly translated into the organismal diversity, i.e. bacterial community diversity per se, because no attempts were made to correct for different copy number of 16S rRNA gene per bacterial genome [Kembel et al., 2012; Stoddard et al., 2015]. Figure 3. Principal components analysis of the data matrix with tree growth and development characteristics as variables and pine climatypes as objects. Active variables: Hei - tree height, Inc - annual increment of tree height, Cro - crone diameter, Dia- trunk diameter, Bra - branching pattern. Supplementary variables: PC1 and PC2 - the 1st and 2nd principal components, extracted from the data matrix with bacterial genera sequence numbers as variables, and rhizosphere soil samples as objects, respectively (i.e. as in Fig.1). Conclusion Here we describe analysis of the bacterial SSU 16S DNA sequence diversity and abundance in the rhizosphere soil of the two different pine species, one of which is a common major forest-forming species (Pinus sibirica) in the Mid-Siberian area where the long-term provenance experiment has been conducted for more than 30 years, the other one (Pinus koraiensis) is a common forest-forming pine species in the Far East. Both pine root systems were dominated by the nucleotide sequence of nitrogen-fixing Bradyrhizobium genus, suggesting its high importance for tree growth. The study showed that pine species influence soil bacterial taxa richness and abundance. The results, providing examples of bacterial taxa with potentially important ecological functions in pine forest ecosystems, extend current knowledge of bacterial community composition in the root-soil interface, which is vitally important for rhizoengineering and rhizoremediation approaches for replanting forest ecosystems.