Hybrid and Environmental E ﬀ ects on Gene Expression in Poplar Clones in Pure and Mixed with Black Locust Stands

: Mixed cropping might be seen as an alternative to monocultures by better protecting biodiversity and improving ecosystem services and resources. In the presented study, we tested the genetic and ecological e ﬀ ects of pure and mixed propagation of di ﬀ erent poplar hybrids planted together with black locust trees. Poplar ( Populus ) hybrids are widely used for bioenergy in monoculture systems due to their rapid and high biomass production. Black locust ( Robinia pseudoacacia L.) is a species with the ability to ﬁx nitrogen and seen as a promising candidate for mixed cultivation. Eight di ﬀ erent poplar hybrids and black locust trees from three provenances planted in two study sites with di ﬀ erent environmental conditions were tested in varying combinations in pure and mixed stands to observe e ﬀ ects of the di ﬀ erent hybrids and genotypes, site conditions and the mixed growing on the performance of poplar and its gene expression. Transcriptome analyses of leaves from four poplar clones selected according to their divergent growth performance were conducted to study di ﬀ erential gene expression that can be an important indicator of di ﬀ erences in growing conditions and success. Di ﬀ erences in gene expression were most pronounced among hybrids and di ﬀ erent genotypes of the same hybrid, followed by the study site inﬂuence, and were least pronounced between mixed and pure stands. The genotypes of the same hybrid were clearly separated from each other. Clear separation between the study sites for all clones was also observed. Only a few genes were di ﬀ erently expressed in pure vs. mixed stand comparisons for each clone, but there were no common genes that were di ﬀ erently expressed in pure vs. mixed stands in all clones. In total, 199 genes showed di ﬀ erential expression between the study sites regardless of poplar clone or type of stands. The analysis suggested that plant genotypes and environmental conditions were more important at the early stage of stand development than pure or mixed cultivation.


Introduction
The growing need for biomass production leads to intensification of land use systems. Most of the economically important species used in agriculture and forestry are grown in monocultures. However, monoculture systems are well-known for their negative consequences for ecosystems resulting in soil degradation and salinization, loss of biodiversity, pollution by fertilizers and herbicides [1][2][3]. In this regard, mixed plant systems based on species with complementary traits could become a sustainable, economically beneficial and ecologically safe alternative or replacement for monocultures [4]. Positive effects of mixed stands for woody species have been demonstrated for Norway spruce, European beech, sessile oak and poplar hybrids [5][6][7][8][9][10].
Poplar hybrids (Populus) are widely used for bioenergy production in monoculture systems due to their rapid growth and high biomass production. Black locust (Robinia pseudoacacia L.) is usually not used for energy production, but is known for its ability to fix nitrogen. Moreover, black locust performs well under drought stress and on poor soils [11,12] and has already been tested in mixed stands with poplar [13][14][15]. Mixing of these two species may potentially improve ecological functions and enhance stand stability.
The fully sequenced and annotated genome of black cottonwood, Populus trichocarpa [16], promoted further transcriptome studies of poplar and other tree species for a better understanding of gene functions of perennial plants in response to their growing environment [17,18]. Since then, many transcriptome studies have been performed with different tree species to investigate their responses to different types of abiotic and biotic stress (e.g., [19][20][21][22][23][24][25]). Importance of genotype by environment interactions for establishment of plantations has been recognized and demonstrated much earlier [4]. Studies conducted with poplar, eucalyptus and black locust species or hybrids [9,11,[26][27][28][29] clearly showed an effect of local conditions on performance and productivity of the genotypes. Many transcriptome studies of plants growing under field conditions were performed for Arabidopsis, rice and maize (see [30] for review), but less for poplar [31][32][33][34].
In the present study, we report the first results for four Populus leaf transcriptomes. Experimental plots were established in two ecologically different sites to observe the influence of the clone, environmental conditions and mixed vs. monoculture cultivation on development of the stands [35]. The following hypotheses were tested in our study: (1) poplar clones respond differently to the planting design; (2) study site conditions influence the performance of stands; (3) the interaction between poplar and black locust trees is weak at the early stage of stand development.

Experimental Design and Sampling
To test for differences in the growth performance of eight commercially used poplar clones in monoculture and mixture with three black locust provenances, poplar cuttings and black locust seedlings were planted at two ecologically different study sites near Göttingen, Germany in April 2014 ( Table 1). The poplar clones and hybrid-clones are well described [36,37] and known for establishing plantations, and were provided by WSD Energieholzbaumschule Björn Diehl (Oberaula, Germany), Staatsdarre Wolfgang (Hanau-Wolfgang, Germany) and P&P Dienstleistungs GmbH Co KG (Eitelborn, Germany). Clonal identity was confirmed by microsatellite analyses. The R3 black locust seedlings (Nagybudmry, Hungary) were purchased under the German legislation of forest reproductive material from P&P Dienstleistungs GmbH Co KG. The two German provenances (R1 and R2) originate from two different provenance regions in Germany and were provided by August Lüdemann (Frankfurt, Germany). One of the study sites, Reinshof, has a young fertile soil with high water storage capacity of Gleyic Fluvisol type (according to the FAO classification), a mean annual temperature of 9.7 • C and an average annual precipitation of 453 mm measured in 2016 by on-site station Adolf Thies GmbH & Co. KG (Goettingen, Germany). The other site, Deppoldshausen, is characterized by the shallow (<60 cm deep) and stony soil with low ability to hold water of Calcaric Leptosol type (according to the FAO classification), a mean annual temperature of 9.0 • C and an average annual precipitation of 457 mm measured in 2016 by on-site station Adolf Thies GmbH & Co. KG (Goettingen, Germany). During the observation period 2015-2016, rainfall and mean annual temperature were lower in Deppoldshausen than in Reinshof (for the weather data see [35]). Four blocks (repetitions) were established at each study site ( Figure S1). The design of the four blocks was exactly the same at each site. Each block consisted of 40 plots (8 rows × 5 columns). The area of each plot was 5 m × 5 m, with 25 trees per plot in total in pure stands and 13 poplar and 12 black locust trees in mixed stands. In each row, one out of the eight poplar clones was planted in a pure plot and in three plots each, mixed with one of the three black locust provenances (Table S1). One of the plots in each row was a pure black locust provenance. Each of the four blocks (repetitions) at one site consisted of the same 40 plots with a random distribution within and between the rows. During the vegetation periods of 2014 and 2015, the growth and survival rates of all poplar clones and black locust provenances were measured [35].
For transcriptome analyses, four clones were selected: two clones "Hybride 275" and "Matrix 11" (P3 and P5, P. maximowiczii × P. trichocarpa) with an intermediate performance in terms of survival and growth, one "Max 1" (P7, P. nigra × P. maximowiczii) with the best performance, and one "Muhle Larsen" (P8, P. trichocarpa) with the worst performance [35]. Only stands with Robinia "HGK 81901" (R1) were selected. The black locust provenances "HKG81901" (R1) and "HKG81902" (R2) showed no significant differences in their performance in terms of mortality and growth [35], thus, they were handled as one provenance throughout the whole experiment. It was observed that bud burst of the P8 clone occurred about two weeks later in comparison to the three other clones. In addition, leaves of clone P8 are covered with oily, strong smelling balsam resin. Fresh young leaves from the middle part of a tree were collected and immediately frozen in the field in liquid nitrogen in May 2016 and stored at −60 • C until RNA extraction. Samples were collected from three clonal trees per plot, in all four repetitions and both study sites resulting in 12 samples per clone in pure and mixed with Robinia "HGK 81901" stands per study site.

RNA Extraction
Total mRNA was extracted from 192 deeply frozen poplar leaf samples (4 clones × 2 sites × 2 types of stands (pure and mixed) × 12 ramets) using the QIAGEN RNeasy Plant Mini Kit (Cat No./ID: 74904; QIAGEN GmbH, Hilden, Germany) and following the manufacturer protocol. Quality of each RNA sample was tested with NanoDrop 2000 (Thermo Fisher Scientific, Waltham, MA, USA) and in 1% agarose gel electrophoresis with 1X TAE as a running buffer. The high quality RNA from five to 12 samples per clone, site and stand representing all four repetitions was pooled to generate 16 mixed equimolar samples that were submitted for RNA sequencing. Some amount of the extracted RNA was stored for further gene expression validation using quantitative real-time reverse transcription PCR (qRT-PCR). In total, 16 mixed samples (4 clones × 2 types of stands × 2 sites) were prepared and sequenced.

RNA Sequencing
Library preparation and sequencing of all samples were done by Chronix Biomedical GmbH (Göttingen, Germany). Additional quality control of RNA was performed with an Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). The cDNA libraries were prepared using the NEBNext Ultra RNA Library Prep Kit for Illumina (New England BioLabs, Frankfurt am Main, Germany). The samples were sequenced using the Illumina NextSeq500 Platform (Illumina, San Diego, CA, USA), which produced read lengths of 75 bp. RNA-seq data have been deposited in the ArrayExpress database at EMBL-EBI under accession number E-MTAB-9041.
The obtained RNA sequencing reads were mapped to the P. trichocarpa transcriptome from Phytozome 12.1 (www.phytozome.net; [38]) using the CLC Genomics Workbench 9.5.4 (CLCbio, Aarhus, Denmark). Trimming for quality and ambiguity of the raw sequence reads were conducted with quality scores (QS) ≥ 30 and default settings of the CLC program. The raw sequence data were processed using fastq [39] with default parameters in R [40]. To check whether the different genetic backgrounds of the clones had an influence on the mapping results and, therefore, on the principal component analysis (PCA) clustering pattern, mapping approaches with differing sensitivity settings (fast and sensitive in Bowtie2 [41], in addition to the default) were tested.

Principal Component Analysis (PCA)
To check for the success of the RNA sequencing and to visualize the grouping of the poplar samples according to clone, type of stand or study site, a PCA was performed. Transcripts with a mean read number below five were not included in the analysis. The PCA was performed using the prcomp function of the stats program and displayed using the ggbiplot program [42] in the R package version 3.4.3 [40].

Identification of Differently Expressed Genes (DEGs) and Sequence Annotation
The analysis of differentially expressed genes (DEGs) among clones, between environments, and pure and mixed stands was performed using the "empirical analysis of DGE" (edgeR) algorithm [43] implemented in the CLC Genomics Workbench program.
Transcripts with at least five mapped reads on average per sample, a false discovery rate (FDR) < 0.1 [44] and a log 2 fold change ≥ 1 were considered to be differently expressed. The poplar gene sequence matches for the DEGs and gene ontology (GO) terms [45] were obtained from Phytozome 12.1 and TAIR (The Arabidopsis Information Resource; http://www.arabidopsis.org) databases [46].

Quantitative Real-Time Reverse Transcription PCR (qRT-PCR)
qRT-PCR was performed to validate differences in gene expression levels. In total, 16 mixed samples were used for the qRT-PCR validation. The cDNA synthesis was performed with 500 ng RNA using the SuperScript III First-Strand Synthesis System for qRT-PCR (Invitrogen, Carlsbad, CA, USA) and Oligo(dT)20 primer (Sigma-Aldrich, Darmstadt, Germany). Three technical replicates for each mixed sample and one housekeeping gene were amplified with a TOptical Gradient 96 Real-Time PCR Thermocycler (Biometra, Analytik Jena, Jena, Germany). Actin 2 was used as a reference housekeeping gene [47].
Two Kunitz trypsin inhibitor genes PtiKPI-A1 (Potri.010G007700) and PtiKPI-A5 (Potri.010G007900) showing strong differences in the expression levels across clones between pure and mixed stands and the study sites were selected for validation (Table S2). The gene specific PCR primers were designed using Primer-BLAST [48] and the corresponding gene sequences from the CLC Genomics Workbench.

RNA Sequencing Data and Transcriptome Annotation
The sequencing of the 16 Table S3). According to the mapping results, reads of P. nigra × P. maximowiczii hybrid P7 clone samples had better mapping to the P. trichocarpa transcriptome than reads of pure P. trichocarpa P8 clone. The mapping was repeated and the results were confirmed.

Principal Component Analysis (PCA)
PCA showed a clear separation among all four clones: the three hybrid clones P3 and P5 (P. maximowiczii × P. trichocarpa) and P7 (P. nigra × P. maximowiczii) and pure black cottonwood species clone P8 (P. trichocarpa) (Figure 1). It was notable that two different genotypes of the hybrid clones P3 and P5 were also clearly separated from each other. A clear separation of the two study sites, Reinshof and Deppoldshausen was observed within all hybrid clusters. The lowest level of separation was observed between pure and mixed stands for each poplar genotype and study site. The same hybrids (P3 and P5) showed a consistent pattern. Clone P7 demonstrated a different pattern in the separation of pure and mixed samples in both study sites: almost no differences were observed between pure and mixed samples in Reinshof, whereas the maximum distance among all other comparisons was observed between pure and mixed samples for clone P7 in Deppoldshausen. Similar or almost identical clustering results were obtained using different sequence matching sensitivity settings for mapping, which indicated that the clustering together of samples from each of the hybrids was not only due to sequence similarity or dissimilarity between the hybrids ( Figure S2).

Differentially Expressed Genes (DEGs) Analyses
Detailed comparisons of DEGs among clones were not performed. This is due to differing sequence similarities between the clones and Populus trichocarpa, which provides the mapping background. Therefore, differences in mapping effectivity would be detected and taken for differences in expression; since these differences will vary in all likelihoods from gene to gene, tests between clones cannot be calculated reliably. For the PCA overview, this was checked by conducting PCAs for different count tables obtained by varying mapping sensitivity settings. On an individual gene level, verification is not possible since full sequence information is not available for all hybrids.
The lack of clustering for mixed vs. pure stands in the PCA was in agreement with the results of the differential gene expression analysis. No common DEGs were identified in the total pure vs. mixed comparison under our experimental conditions and limited numbers of replicates. Nevertheless, this reflects low effects of the heterospecific neighbor (black locust) in the poplar stands at this stage of plantation development and high genetic variability between clones. Single clone pure vs. mixed comparisons for each study site separately could not be performed due to the lack of replicates within each genotype for sufficient statistical power to identify clone-specific interactions. By combining the samples from the two study sites, we could perform the comparison of pure vs. mixed stands for each single clone. These revealed genes that are influenced by the type of stand independently from the study site conditions. In total, 30 genes showed differential expression between pure and mixed samples ( Table 2). Only two genes were common in the data set across clones-Potri.017G025900 in P7 and P8 with opposite expression behavior, and Potri.001G441400 in P3 and P7. In each comparison, stress and/or defense response genes were detected as DEGs and showed lower expression levels in the mixed stands.
According to the PCA, the second strongest clustering after the differences between hybrids was the clustering for study sites. To evaluate the influence of a study site on clone performance we compared the complete sample sets of Deppoldshausen vs. Reinshof. In total, 199 genes showed differential expression between the two study sites in combined comparisons regardless of clone or type of stands (pure vs. mixed) (Table S4). Among them, 119 genes showed higher expression levels in Deppoldshausen and 80 in Reinshof.
The study site comparisons were also performed for each clone separately. The highest number of DEGs between the study sites was observed for clone P7 with total of 1706 genes, and similar numbers were observed for other clones (Table 3, Tables S5 and S6).  The clone-specific differences in DEG numbers between Deppoldshausen and Reinshof were much stronger than in the total study site comparison. For clone P3, 1207 DEGs were found between the study sites, of which 685 genes showed higher expression levels in Reinhof and 522 in Deppoldshausen.
For clones P5, P7 and P8, an opposite pattern was observed: about two-thirds of all DEGs between the study sites showed higher expression levels in Deppoldshausen than in Reinshof. A low number of genes in the total study site comparison and much higher gene numbers in the clone-specific comparisons showed that individual genotype has higher influence on the transcriptome response of plants than environmental conditions.

Quantitative Real-Time Reverse Transcription PCR (qRT-PCR)
Two genes from the Kunitz trypsin inhibitor family PtiKPI-A1 (Potri.010G007700) and PtiKPI-A5 (Potri.010G007900) were selected for qRT-PCR of mixed samples to validate the RNA-seq data. The qRT-PCR confirmed expression profiles for all clones in the mixed RNA samples (Table S7).

Discussion
The aim of this study was to detect environmental effects on gene expression among poplar clones with different field performance and possible signs of interactions between poplar and black locust plants. The plots were established in April 2014. At the time of sampling in May 2016, plants had two years to adapt to the site conditions and develop sufficient root systems.
The use of relatively short 75 bp long single-end reads in our study is well-justified by simulation and real-study experiments, which showed that 50 bp reads were sufficient for good quality mapping and DEG detection and that results based on 75 bp long reads were similar to those based on 100 bp long reads [52,53]. Moreover, 50 bp long reads and poplar transcriptome annotation from Phytozome were successfully used in poplar stress response studies [54,55].
It is interesting that reads of P. nigra × P. maximowiczii hybrid P7 clone samples had better mapping to the P. trichocarpa transcriptome than reads of pure P. trichocarpa P8 clone, but the difference was not really great and was likely occasional due to high genetic similarity between these poplar species, especially P. maximowiczii and P. trichocarpa, which belong to the same Populus section, namely Tacamahaca [56].
The PCA analysis showed clear differences between the transcriptomes of hybrids, genotypes, study sites and type of stands (pure vs. mixed stands). All three hybrids were separated from each other based on their species background. The two clones of the P. maximowiczii × P. trichocarpa hybrid (P3 and P5) clustered together close to clone P7 (P. nigra × P. maximowiczii), while clone P8 (P. trichocarpa) was completely separated from the three other clones. Different poplar species and their hybrids have been successfully distinguished based on various phenotypic traits [57] and molecular markers, such as the chloroplast trnT-trnF region and nuclear rDNA sequencing [58], amplified fragment length polymorphisms (AFLPs) [59], microsatellites [60] and single nucleotide polymorphisms (SNPs) [61]. Our results demonstrated that a very clear taxonomic separation is also possible using transcriptome data. Two study sites, Reinshof and Deppoldshausen, were also clearly separated from each other for all hybrids. The weakest differentiation was observed between pure and mixed stands. Clones P3, P5 and P8 showed a similar pattern in the study site and type of stands clustering, but for clone P7, clear differences in separation of pure and mixed samples according to the study site were detected. The samples from Reinshof, as the site with favorable growing conditions, were almost overlapping. At the same time in Deppoldshausen, as the marginal growing site, very clear differences between pure and mixed samples were observed. This suggests that clone P7, despite of its best growing performance, is more sensitive to the growing conditions in comparison to other clones. In summary, the PCA analysis revealed that environmental conditions and genotypes had a more pronounced influence than the type of stand (pure vs. mixed) on differentiation between poplar clones at the transcriptome level, at least at this early stage of stand development.
The differential gene expression analysis based on the RNA-seq data was in agreement with the results of the PCA analysis. Total comparison between pure and mixed samples did not reveal any DEGs. This finding matched very well with basal area measurements and mortality data of poplar trees [35]. No significant differences in basal area between pure and mixed stands were observed. Likewise, total mortality data for poplar also showed no significant differences between mixed and pure stands in 2014-2015. However, significantly higher mortality in the mixture than in pure stands was detected to the end of 2016, which indicated increasing competition pressure of black locust. Only 30 genes showed significant differences in expression between pure and mixed samples for all four clones regardless of the study site ( Table 2). DEGs involved in stress and defense response were observed for each clone. All of them showed lower expression levels in mixed samples. Previous studies showed that there seem to be universal stress response genes that always respond to any kind of stress [62,63]. Some studies demonstrated stress reduction in mixed stands during biotic [64] or abiotic [65] stress. Therefore, we hypothesize that the lower expression levels of stress and defense response genes in our mixed stands indicates that the slightly increased biodiversity (two species instead of one) already reduces stress to the plants. Mortality, growth and RNA-seq data together showed weak interaction between poplar and black locust, as a nitrogen-fixing species, in the first years of stand development (from 2014 until the start of 2016). Rooted seedlings of black locust were in advantage in growth over poplar trees that were planted as unrooted cuttings and had to develop a sufficient root system. Fertile soil conditions also promoted the growth of black locust and delayed development of interactions between both species, even in Deppoldshausen with its poorer than in Reinshof soil conditions. The same effect has been observed in studies with other N-fixing species [66,67].
Reinshof and Deppoldshausen were clearly separated by 199 DEGs regardless of the type of clone or stand (Table S4), with more genes showing higher expression levels in Deppoldshausen. The clone comparisons for gene responses to contrasting study sites were in agreement with the results of the PCA. Deppoldhausen and Reinshof were clearly separated in their effects on gene expression in all clones and also demonstrated clonal differences (Table 3). About two-thirds of all DEGs between the study sites showed higher expression levels in Deppoldshausen than in Reinshof. Clone P7 had the highest number of DEGs, most of them showing higher expression levels in Deppoldshausen. Max 1 hybrid (P7) is known for its very good growing performance and high survival rates which might explain higher DEG values [36,37]. The marginal growing site Deppoldshausen provoked a stronger genetic response than the fertile growing site Reinshof. The separation between clones and sites has also been observed for DEGs in wood [32].
After two years of the stands' development, weak differences were still observed between pure and mixed stands and higher importance of growing sites for tree performance. The effects of species genetic background were the strongest. Based on this study, the clone Max 1 (P7) was selected as the most promising and interesting candidate for further transcriptome studies. According to our results, it combines the best growing performance with sensitivity to growing site conditions.

Conclusions
In this study, we report the first results on transcriptome analyses in poplar hybrids cultivated under different field conditions in pure stands and mixed with black locust stands. The results of the transcriptome analysis showed the clearest and the strongest differences between hybrids and genotypes of the same hybrid, followed by the study sites; the smallest differences were between types of stand. No DEGs were detected in the total comparison of pure and mixed stands as the result of high genetic clone specific variability of poplar clones. The small number of DEGs revealed for each clone in pure vs. mixed comparisons confirmed low interaction between poplar and black locust at that early stage of development, because it is likely that trees were too young, and their root systems had not formed any interactions in soil yet. There were 199 DEGs between the two study sites in total regardless of poplar clone or type of stands, which showed the importance of study site conditions for tree growth. The number of DEGs differed between clones depending on the study site. In general, the results indicate a strong influence of local environmental conditions and genotypes of clones on the performance of trees in pure and mixed stands at this early stage of the project. Further studies should be carried out with the samples based on preferably single clones to exclude genotype effects, but with a sufficient number of biological replicates for more detailed expression analysis.
Supplementary Materials: The following are available online at http://www.mdpi.com/1999-4907/11/10/1075/s1, Table S1: Plot design with eight poplar clones and three black locust provenances, Table S2: Kunitz trypsin inhibitor genes PtiKPI-A1 and PtiKPI-A5 and the corresponding PCR primers used for qRT-PCR validation, Table S3: The quality and quantity of sequenced reads for all samples before and after trimming, mapping and coverage values for all poplar samples, Table S4: List of 199 DEGs in the both study sites observed in all poplar clones and pure and mixed stands, Table S5: List of DEGs in Deppoldshausen for each clone in comparison between the study sites, Table S6: List of DEGs in Reinshof for each clone in comparison between the study sites, Table S7: Differences in expression based on RNA-seq and qRT-PCR of two genes from the Kunitz trypsin inhibitor family in mixed RNA samples between pure and mixed stands of different poplar clones in two study sites: Reinshof and Deppoldshausen, Figure  Funding: This research was funded by the German Federal Ministry of Education and Research (BMBF; FKZ 031A351C) within the framework of the IMPAC 3 project ("Novel genotypes for mixed cropping allow for improved sustainable land use across arable land, grassland and woodland").