Next Article in Journal
Stem Decay in Live Trees: Heartwood Hollows and Termites in Five Timber Species in Eastern Amazonia
Previous Article in Journal
Comparison of Two Sampling Methods to Estimate the Abundance of Lucanus cervus with Application of n-Mixture Models
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Forests
Volume 11
Issue 10
10.3390/f11101086
forests-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Genetic Consequences of Hybridization in Relict Isolated Trees Pinus sylvestris and the Pinus mugo Complex

by
Karolina Sobierajska
1,†,
Witold Wachowiak
2,†,
Julia Zaborowska
2,
Bartosz Łabiszak
2,*,†,
Błażej Wójkiewicz
1,
Maciej Sękiewicz
1,
Anna K. Jasińska
1,
Katarzyna Sękiewicz
1,
Krystyna Boratyńska
1,
Katarzyna Marcysiak
3 and
Adam Boratyński
1
1
Institute of Dendrology, Polish Academy of Sciences, 62-035 Kórnik, Poland
2
Institute of Environmental Biology, Adam Mickiewicz University in Poznań, 61-614 Poznań, Poland
3
Department of Botany, Kazimierz Wielki University, 85-064 Bydgoszcz, Poland
*
Author to whom correspondence should be addressed.
†
These authors contributed equally to the work.
Forests 2020, 11(10), 1086; https://doi.org/10.3390/f11101086
Submission received: 17 September 2020 / Revised: 5 October 2020 / Accepted: 6 October 2020 / Published: 12 October 2020
Browse Figures
Versions Notes

Abstract

:
Scots pine (Pinus sylvestris L.) and the taxa from the P. mugo complex can hybridize in the contact zones and produce fertile hybrids. A unique example of an early Holocene relict population of P. sylvestris and P. uliginosa (a taxon from the P. mugo complex) growing on the tops of Jurassic sandstone rocks is located in Błędne Skały (Sudetes). Phenotypically, there are trees resembling P. sylvestris, P. uliginosa and intermediate forms between them. We expected that some of P. sylvestris and/or P. uliginosa-like trees could be in fact cryptic hybrids resembling one of the parental phenotypes. To address this question, we examined randomly sampled individuals, using a set of plastid (cpDNA), nuclear (nDNA) and mitochondrial (mtDNA) markers as well as biometric characteristics of needles and cones. The results were compared to the same measurements of allopatric reference populations of the P. sylvestris and the P. mugo complex (Pinus mugo s.s, P. uncinata and P. uliginosa). We detected cpDNA barcodes of the P. mugo complex in most individuals with the P. sylvestris phenotype, while we did not detect cpDNA diagnostic of P. sylvestris within P. uliginosa-like trees. These results indicate the presence of cryptic hybrids of the P. sylvestris phenotype. We found only three typical P. sylvestris individuals that were clustered with the species reference populations based on needle and cone characteristics. Most trees showed intermediate characteristics between P. sylvestris and P. uliginosa-like trees, indicating intensive and probably long-lasting hybridization of the taxa at this area and subsequent gene erosion of parental species.

1. Introduction

Hybridization and introgression are important evolutionary factors that increase species variation and may lead to speciation [1,2,3,4,5,6]. These two processes concern mostly related species, which do not express isolation mechanisms such as phenological barriers [7] and post-pollination barriers to fertilization [8]. The interbreeding between different species of Pinaceae was found to be quite frequent [9,10,11,12]. Considering their long generation time and lifespan, the relatively recent divergence between dwarf mountain pine (P. mugo Turra) and Scots pine (P. sylvestris L.) [13] resulted in their similar genetic variation at neutral loci [14]. The hybrids between the taxa of the P. mugo complex and P. sylvestris in the wild state have been traditionally reported as morphologically intermediate individuals ([15,16] and the literature cited herein). The biometric studies were conducted to verify the presence of hybrid individuals or even populations [17,18,19,20,21,22,23,24,25]. However, the unequivocal identification of hybrids is difficult based on morphological and anatomical characteristics alone.
During the last few decades, isoenzymatic and molecular analyses were used to estimate the potential hybridization between the taxa of the P. mugo complex and P. sylvestris [26,27,28,29,30,31,32,33]. The possibility of hybridization was also investigated by artificial crossing experiments ([10,34,35] and the literature cited herein). The results indicated the higher success of pollination, when taxa of the P. mugo complex were the pollen donors, and possible back-crossing of hybrids with both parental species. However, a more complete picture of interspecific gene flow within natural populations of hybridizing taxa could be obtained by combining biometric, biochemical and molecular methods [19,30,36,37,38].
To date, studies on hybridization between the P. mugo complex and P. sylvestris have been conducted in different types of the species contact zones including the small and isolated populations of one taxon surrounded by extensive, massive populations of the other [10,25,26,30] or between relatively large populations of both species [31]. However, to our knowledge, hybridization has not been investigated in small, relict and isolated populations of the taxa. Such stands are rare but are especially valuable, as relict populations may conserve products of hybridization and introgression repeated over many generations without or with only very restricted gene inflow from outside. A population of this type is located in Błędne Skały in the Stołowe Mountains (Sudetes). Pines grow there on the tops of Jurassic sandstone rocks, in hardly accessible places, which have not been exploited for wood and are unsuitable for pasturing [39,40]. This mixed stand of pines originated most probably from the Boreal period of the Holocene. The pollen of P. sylvestris type, which may also include very similar and indistinguishable pollen of the P. mugo complex [41], appeared during that time in high frequencies in the Sudetes [42,43], and was detected also on peat-bogs nearest to the Stołowe Mountains [44,45]. The early Holocene expansion of broad-leaved trees and Norway spruce (Picea abies L.), which formed dense and shady forests, reduced occurrence of light demanding P. sylvestris and P. mugo to the sites where conditions were inaccessible for other forest tree species such as Corylus, Ulmus, Alnus and Picea [44,45,46]. Consequently, pines could survive in the Sudetes only on the rocks at the tops of the mountains and/or the rocky edges of their plateau. The best examples of such sites in the Stołowe Mountains are the sandstone labyrinths on the tops of Szczeliniec Wielki and Skalniak mounts. Regarding Błędne Skały location, there are currently growing individuals phenotypically resembling P. sylvestris, P. mugo s.s. and P. uliginosa G. E. Neumann [47]. The pollen record of P. sylvestris type proves the persistent character of pine population in the Stołowe Mountais during the entirety of the Holocene [44,45,48], and even reflects the temporal expansion of pines after fires in the High Middle Ages and Early Modern period (see Table 2 and Figure 5 in [48]).
Pines have high pollen production and dispersal ability [49,50,51], however, most pollen grains fall in a distance up to a few hundred meters from the source tree [52]. Therefore, considering the species biology, success of pollination and the recognized mating system [52,53], the inflow of pollen from outside Błędne Skały population was most likely very limited from the boreal time in the case of rare P. uliginosa (the species from the P. mugo complex) and at least significantly reduced in the case of much more common P. sylvestris. Thus, currently the observed mixed population of both taxa, which consist of individuals phenotypically resembling typical species and intermediate forms, was formed and has persisted locally over the last 7000–8000 years. In this context, the population in Błędne Skały is a unique and very interesting model for hybridization and gene flow studies.
In the presented study, we investigated the hybridization processes in the pine population in Błędne Skały to verify (i) if there are any hybrids among P. sylvestris-like and/or P. uliginosa-like trees which resemble one of the parental trees; (ii) if there are any patterns of successful hybridization as compared to some earlier works (iii) to what extent the level of genetic variation of this small and isolated population is comparable to variation in much bigger reference populations of pure species. To address those questions, we sampled individuals across the phenotypic forms present in this population and performed genetic and biometric analyses to verify taxonomic relationships between the samples and to compare morphology of trees from the mixed population to the reference populations.

2. Materials and Methods

2.1. Material

The mixed pine population in Błędne Skały (BS) occurs on the upper parts of sandstone rocks and is located at the summit of the north-western part of the Skalniak ridge in the Stołowe Mountains, at an elevation of 850 m a.s.l. It consists of dispersed pine individuals of P. sylvestris and P. uliginosa phenotype, with a few individuals resembling P. mugo s.s. and a number of individuals of intermediate phenotype. The pine population is surrounded by extensive Picea abies forests, which isolate it from direct contact with other pine populations.
Material for the study was taken randomly from 64 individuals in an area of about 2 ha. The sampled individuals were about 20 m distance from each other or growing on separate rocks. During sampling, part of the individuals was preliminary identified in the field based on growth form, bark color on the upper part of stem, needle color and shape and setting angle of the one-year-old conelet [15,16]. The sampling included 15 individuals of P. sylvestris-like phenotype (BS_S), 17 individuals resembling the taxa from the P. mugo complex namely P. uliginosa (BS_M), 7 individuals of intermediate character between P. sylvestris and P. uliginosa classified as hybrids (BS_H) and another 25 unclassified individuals (BS_N). As the sampled trees could potentially represent a genetic background of different pine species, we included a set of reference populations taken from across current natural ranges of P. sylvestris and the taxa of the P. mugo complex (P. mugo s.s., P. uliginosa and P. uncinata Ramond (Table 1; Figure 1). These populations have been used in previous morphometric [54,55,56,57,58,59,60,61] and genetic [30,62,63,64,65] studies. Particular analyses were conducted using all taxa, but not always the same reference populations.

2.2. DNA Extraction, Amplification and Marker Genotyping

Genomic DNA was isolated from dry needles, using a DNeasy Plant Mini Kit (Qiagen, Hilden, Germany). DNA quality was evaluated using a BioPhotometer (Eppendorf, Hamburg, Germany). The samples were analyzed with a set of organellar and nuclear markers that show different modes of inheritance. Polymorphisms at 14 mitochondrial regions, inherited in pines in the maternal line and distributed at short geographical distances by seeds, were investigated including 13 markers [66], and the insertion–deletion (IN–DEL) locus nad1 intron B/C [67]. The standardized 15 µL PCR reaction mix contained 1 × PCR buffer, 1X bovine serum albumin (BSA), 1.5 μM MgCl2, 40 µM dNTP mix, 0.2 µM of both forward and reverse primers (Genomed, Warsaw, Poland), 0.15 U Taq DNA polymerase (Novazym, Poznań, Poland) and 15–20 ng of DNA template. Reaction conditions were unified with 3 min of initial denaturation at 94 °C, 35 cycles with 30 s denaturation at 94 °C, 30 s annealing at 57°C and 1 min 20 s elongation at 72 °C, with the final extension for 5 min at 72 °C. Variants were genotyped by polymerase chain reaction restriction fragment length polymorphism (PCR-RFLP) and digested in 10 µL reaction volume, containing 3 µL of the PCR product, 2X buffer and 2.5 U of the relevant restriction enzyme [68]. Samples were incubated for at least 10 h, following the manufacturer´s recommendations for respective enzymes (Thermo Fisher Scientific, Waltham, MA, USA), and electrophoretically separated on 1.5%–2% agarose gel, stained with GelRed in UV. The insertion–deletion polymorphisms in intron B/C of the nad1 locus was amplified and scored according to the methods described in Soranzo et al. [67]. The region contains shorter variant a that appears to be mostly fixed in P. mugo s.s. populations and it segregates with different frequency across the geographical distribution of P. sylvestris. The polymorphisms at nad1 B/C region was scored by 2% agarose gel electrophoresis of the corresponding PCR-amplified fragment.
Variation at the plastid genome, which is paternally inherited in these taxa and transmitted by pollen, was studied at a set of 12 microsatellite markers. Two multiplex reactions of the chloroplast simple sequence repeats (cpSSR) loci were run, including regions PCP1289, PCP26106, PCP30277, PCP36567, PCP41131, PCP45071, PCP87314, PCP102652 (multiplex 1) and Pt15169, Pt26081, Pt30204, Pt71936 (multiplex 2) [69]. Additionally, a species-diagnostic plastid (cpDNA) intergenic trnF-trnL region that discriminates P. sylvestris from the taxa of the P. mugo complex was used to determine the origin of the plastid genome in samples from the focal population in Błędne Skały. This PCR-RFLP marker has previously been developed based on a single-nucleotide polymorphism that leads to an undigested PCR product for P. sylvestris and a digested PCR product for P. mugo when treated with DraI enzyme [70]. The PCR was run in a total volume of 15 μL, containing 15 ng of template DNA, 10 μM of each of dNTP, 0.2 μM each of forward and reverse primer, 0.15 U Taq DNA polymerase, 1X BSA, 1.5 μM of MgCl2 and 1X PCR buffer (Novazym, Poznań, Poland). Standard amplification conditions were applied, including initial denaturation at 94 °C for 3 min, followed by 35 cycles with 30 s denaturation at 94 °C, 30 s annealing at 60 °C and 1 min 30 s extension at 72 °C, with a final 5 min extension at 72 °C.
Variation at the biparentally inherited nuclear genome was investigated using a set of nine simple sequence repeats (nSSR), including psyl2, psyl25, psyl36, psyl42, psyl44, psyl57 and ptTX3025, ptTX4001, ptTX4011, following two multiplex reactions [71]. The fluorescently labelled PCR products, obtained for the studied cpDNA and nuclear (nDNA) microsatellite regions, along with a size standard GeneScan 500 LIZ (Applied Biosystems, Foster City, CA, USA), were separated on a capillary sequencer ABI 3130 (Thermo Fisher Scientific, Waltham, MA, USA). The SSR allele size was determined using the GeneMapper software ver. 4.0 (Applied Biosystems, Foster City, CA, USA).

2.3. Diversity and Differentiation Measures

Alleles identified within each mitochondrial (mtDNA) marker were concatenated, and connections among the obtained haplotypes were analyzed by a Median Joining network constructed in PopART 1.7 [72,73]. The number of segregating sites (S), number of private variants (SP) and the mean within group distance (d) were calculated in DnaSP V6 [74] and MEGA 7 [75] for groups of samples defined as BS and for the reference populations. Similarly, the numbers of haplotypes (Nh) and private haplotypes (Np) present in populations and species were calculated. To better assess the scale of genetic variation at mtDNA segregating in the taxa, the haplotype frequencies and population gene diversities (Hd) were estimated in DnaSP. The relationships among samples were evaluated based on the mean distances between populations (dxy), using the Unweighted Pair Group Method with Arithmetic Mean (UPGMA) as implemented in MEGA. The distribution of genetic variation among populations, and among genetic structures detected in the population groups were checked by Fst and Gst calculated in DnaSP. The differences between populations were verified by hierarchical analysis of molecular variance (AMOVA), using Arlequin 3.5.22 [76].
To describe the variation of cpDNA SSR loci, for each defined sample group we calculated mean number of haplotypes (Ah), mean effective number of haplotypes (Ae), number of private haplotypes (AP), haplotype richness (AR), haplotype diversity (Hd) and mean genetic distance of individuals within populations (D2sh) using GenAlEx 6 [77]. Variation at the diagnostic trnF-trnL cpDNA region was scored as a presence or absence of a DraI restriction site.
Genetic variation at nSSR of the mixed population and the reference populations was analyzed based on mean number of alleles per locus (Na), mean effective number of alleles per locus (Ne), number of private alleles (Ap), mean observed heterozygosity (Ho), mean expected heterozygosity (He), mean fixation index (F). The analysis of molecular variance (AMOVA) on the Rst distances, as implemented in Arlequin 3.5.22, was used to evaluate the differentiation among defined groups of populations from BS contact zone and reference populations at nSSR loci. Afterward, the principal coordinates analysis (PCoA) on Nei’s genetic distance was applied to visualize the genetic structure of the populations and BS group. The Bayesian clustering method implemented in STRUCTURE 2.3.4 [78,79,80] was used to assign individuals and populations into genetically distinct groups without a priori assumption of their origin. The STRUCTURE software uses the algorithm of Markov chain Monte Carlo (MCMC) to assign particular individuals to a number of genetic clusters (K) whose members share similar patterns of variation. It does not consider samples’ origins but assumes that each cluster is in optimal Hardy–Weinberg and linkage equilibria. The parameter sets, that assumed correlated allele frequencies and the used admixture model, together allowed attribution of mixed recent ancestry of individuals and assignment of the proportion of inferred clusters’ origin within each individual genome. The “recessive alleles” option was used, with 20 runs for each K, from K = 1 to 10, and with burn-in lengths of 500,000 and 100,000 iterations. The data probability distributions (Ln P [D]) and the ΔK values [81] were calculated using STRUCTURE HARVESTER Web v0.6.94 application and were used to determine optimal K for given dataset [82].

2.4. Biometrics

We characterized two-year-old needles using a set of traits (Table S5), following the procedures described by Boratyńska et al. [56] and Boratyńska and Boratyński [23]. Every individual was characterized using the average measurements of five needles, each from a different dwarf shoot (i.e., the morphological structure in Conifers from which needles in fascicles of 2–5 grow together). The anatomical characteristics were measured and/or evaluated on the preparations from the central part of the needle [21,23,55,56,83].
Cones were characterized using 11 measured characteristics and five proportions (Table S6) describing cone and cone scale shape [17,54,84]. Measurement procedures were adopted from Marcysiak [54] and Marcysiak and Boratyński [57]. Length, width, diameter and circumference measurements were taken from closed cones after soaking, while the number of cone scales was determined for open cones after drying. The scale characteristics were measured on one scale taken from the convex side of each cone at the maximum diameter. All measurements were done manually, using an electronic caliper.

2.5. Statistical Analyses of Biometric Data

Prior to multivariate comparisons, the data distribution and homoscedasticity of the variances were inspected, applying the Shapiro–Wilk’s and Brown–Forsythe tests, respectively [85,86]. The arcsine transformation of percentage data was applied prior to statistical analyses. The data were standardized using the STATISTICA 9 software (StatSoftPolska, Kraków, Poland) to avoid a possible influence of their different types on the results of the multivariate analyses.
We determined the minimum and maximum values of the characteristics, arithmetic means, and variation coefficients. These parameters were calculated and analyzed for every individual in the BS population, and for reference populations of P. sylvestris, P. mugo, P. uliginosa and P. uncinata. Principal component analysis (PCA) was used to detect grouping of samples from the BS mixed population and their position among referenced populations of P. sylvestris, P. mugo and P. uncinata. Subsequently, discrimination analysis was conducted on these characteristics, except for the ones which were strongly biased (SRC, SRF, VB and RC – see Tables S5 and S6 for traits description) [86,87].
To assess the significance of differences for specific measured traits between groups preliminary determined in the BS population during sampling, Tukey’s HSD test for traits with normal distribution and the Kruskal–Wallis test for traits with biased distribution were applied [86,88], using STATISTICA. The taxonomic position of each individual from the BS mixed population, assessed via macro-morphological traits during sampling, was verified using characteristics of the needles and cones. Finally, the results from the biometric study were compared with the results of the genetic analyses to distinguish individuals of P. sylvestris from the taxa from the P. mugo complex and their hybrids.

3. Results

3.1. Genetic Diversity and Differentiation

In total, 39 mitotypes were identified across the four pine species (Figure 2, Table S1). Species-specific combination of haplotypes was observed only in P. uncinata, which was distinct when compared to the other species (Figure 2, Figure S1). Different frequencies of mitotypes were observed in individual reference populations. The highest frequency of shared haplotypes was found in individuals from the mixed BS population (Table S2). In the BS pine mixed population, we observed mitotypes of P. sylvestris, P. mugo s.s. and P. uliginosa reference populations, but not of P. uncinata. Unique mtDNA haplotypes of pines from BS population were found for five individuals from hybrid or unclassified sample groups. The remaining pines shared numerous mitotypes present in individuals from BS and reference populations. The most frequent mitotypes at the BS population were found at high frequency in the P. uliginosa population from Węgliniec Nature Reserve (Table S1). Some unique mitotypes were found at low frequency in every examined species and BS population. The samples including four reference species and a BS population formed each separate group in the UPGMA tree based on Nei’s genetic distance at mtDNA regions (Figure S2).
Only three individuals from the BS population had diagnostic a cpDNA marker characteristic for P. sylvestris. Those three individuals were indicated as pure P. sylvestris based on STRUCTURE analysis of cpSSR loci. The remaining individuals from BS had plastid DNA specific to the P. mugo complex (i.e., P. uliginosa), independently of the phenotype.
Individuals from BS showed reduced genetic variation at cpSSRs, with a considerably lower number of haplotypes and lower haplotype diversity as compared to the reference species (Table 2A) and individual populations (Table S3). The individuals unclassified taxonomically and those defined as P. sylvestris from the BS mixed pine population showed the highest mean genetic distance between individuals within populations (D2sh range 32.1–51.9) as compared to the other groups in this area (7.3–9.1) and the reference populations (4.4–13.6). Haplotype and gene diversity at mtDNA markers were also lowest in the BS population (Table 2B). However, in contrast to organelle markers, variation at biparentally inherited nSSR loci was similar in the BS population as in the reference samples (Table 2C, Table S4). There was no evidence of reduced heterozygosity among samples in the BS population (Table 2C).
According to the PCoA based on genetic distance at nSSR loci, all defined groups of individuals from BS showed close genetic similarity to each other, to P. uliginosa from the Węgliniec Nature Reserve and to P. uncinata samples (Figure 3). This group of samples was distinct from the reference populations of P. mugo s.s. and P. sylvestris. Among two outlier populations detected appeared one P. uliginosa from the Large Batorów Peatland, with high genetic similarity to group of three P. mugo s.s. populations and the fourth P. mugo stand localized in Sudetes (M8), which was distinct as compared to all other samples analyzed (Figure 3).
Structure analysis at nSSR loci indicated the presence of five genetically distinct groups corresponding approximately to the four reference species and one sample group from BS (Figure 4, Figure S3). However, the pattern was not homogenous, and the samples showed varying levels of admixture among different clusters. Three genetic clusters, including P. sylvestris, P. mugo s.s. and a mixture of samples from populations of other taxa, were found at cpSSR loci (Figure S4).

3.2. Phenotypic Variability and Differentiation

The level of the intraspecific morphological variation differed among the taxa and populations from BS in terms of needle (Table S5) and cone characteristics (Table S6). The highest variability in most of the needle characteristics showed P. uliginosa and the population from BS; reference taxa variation was at a lower level (Table S5). Regarding the cone characteristics, the most variable ones were P. sylvestris and P. uliginosa (Table S6).
Most of the analyzed needle and cone characteristics differed significantly (p ≤ 0.01) between samples from BS and reference populations. The BS population significantly differed in terms of needle characteristics from every reference species, while in terms of cone characteristics it differed significantly from P. uncinata and P. mugo (p ≤ 0.01). The groups of individuals described in the putative hybrid population during material collection based on macro-morphology as P. sylvestris-like (BS_S), P. uliginosa-like (BS_M), hypothetical hybrids intermediate between these two taxa (BS_H) and individuals not determined (BS_N), appeared close to each other in terms of needle and cone characteristics. Each of these four groups did not reveal significant pairwise differences neither in terms of needle nor in terms of cone characteristics. The individuals determined during sampling as P. uliginosa-like (BS_M) in terms of needle characteristics revealed a relatively low number of traits differentiating them from reference populations of that species, but a relatively high number of characteristics differentiating them from P. mugo s.s. and P. uncinata. Twice as many characteristics differentiated the group BS_M from P. sylvestris populations (Table S7). However, there was no general rule in terms of cone characteristics (Table S8). Surprisingly, there were no differences between each group from BS and the P. uliginosa population from the Zieleniec Nature Reserve (UL2).
The centroids of four groups from BS were close to each other on the scatterplots of the PCA and fell between reference populations of P. sylvestris and the taxa representing the P. mugo complex.
The group assigned to P. sylvestris during sampling appeared only slightly similar to these species reference populations, and only in terms of needle characteristics (Figure 5). The needle characteristics that mostly contributed to group differentiations were NRC, SC, SF, TN/WN, VBF and VBL, while the most important cone characteristics were CD, CL, CDM, DAU, TA, CVX and WA (Figure S5, code descriptions in Tables S5 and S6).
The needle characteristics on the dispersion diagram between the two first discrimination variables, responsible for 80% of the total variation, indicate the intermediate position of BS samples between P. uncinata and individuals of P. uliginosa and P. sylvestris. Some of BS individuals intermixed with individuals of P. uncinata, P. uliginosa and P. mugo, and entered the 95% confidence intervals of these species. Interestingly, BS individuals did not intermix with P. sylvetris specimens (Figure S6A). In the cone characteristics, pines from BS partly intermixed with reference P. uncinata, P. uliginosa and P. sylvestris, but did not enter P. mugo s.s. at a 95% confidence interval (Figure S6B).
The dispersion of BS samples described by the two first PCA variables also indicates their intermediate character among the compared taxa. In the needle characteristics, a great part of BS individuals, determined during sampling as P. uliginosa-like, was located between P. mugo s.s., P. uliginosa and P. uncinata, while the remaining ones were located between P. uliginosa and P. sylvestris. Individuals determined as BS_S were located between reference P. sylvestris and taxa of the P. mugo complex, while those determined as hybrids were dispersed among all the others (Figure 6A). In the cone characteristics, BS_M individuals mostly revealed close connection to reference P. uliginosa, P. uncinata and P. mugo s.s., although some of them were located between P. sylvestris and taxa from the P. mugo complex, with high affinity to P. uncinata.
In PCA analysis, trees from BS assigned to P. sylvestris were located between reference P. sylvestris and P. uncinata populations, while the putative hybrids were more dispersed (Figure 6B). The three individuals identified as pure P. sylvestris in genetic analysis were close to this species´ reference populations on the PCA, both based on needle and cone characteristics (Figure 6) and in the discrimination analysis (Figure S6).

4. Discussion

Due to the relatively recent speciation history and the possible gene exchange during speciation [13], the taxa of the P. mugo complex and P. sylvestris show high genetic similarity at nuclear and organellar genomes [62,63,64,71]. So far, species diagnostic markers were found at cpDNA discriminating P. sylvestris from the taxa of the P. mugo complex [33,35,70]. In advances on earlier studies, we provide another diagnostic marker at mtDNA, which discriminates between P. uncinata and the remaining species from the P. mugo complex. However, the unique P. uncinata mitotype was not observed in the studied BS population, indicating that P. uncinata did not contribute, at least as a seed donor, to the development of the group of pines at the study area. Pines from the BS stand had most of the mtDNA mitotypes observed at high frequency in the small, relict population of P. uliginosa surrounded with extensive P. sylvestris stands in the WęgliniecNature Reserve. However, the discriminating power of the mtDNA is still too low to be conclusive about any further genetic relationships between the samples and between P. mugo s.s. and P. uliginosa.
The indication of a hybrid origin of the individuals comes from the analysis of the paternally inherited plastid and the biparentally inherited nuclear genome [30,62,64,71]. Interestingly, only three individuals carrying the cpDNA genome typical for P. sylvestris were identified among the 64 individuals sampled: one of them was assigned based on the macro-morphological phenotype to P. sylvestris (BS_S), while the other two remained undetermined during sampling (BS_N). All remaining individuals analyzed from BS possessed plastid DNA characteristic of the P. mugo complex, independently of phenotype.
The plastid genome derived from the P. mugo complex observed in the 14 out of 15 trees classified phenotypically as P. sylvestris means that those samples represent cryptic hybrids, which resemble pure species and cannot be distinguished in the field via macro-morphological characteristics. Furthermore, despite the vast diversity of phenotypic forms sampled from Błędne Skały and preliminarily classified into three distinct phenotypic groups, they all showed a uniform chloroplast background and high genetic similarity to each other at nuclear loci. For instance, they had the lowest genetic distance and high similarity at the allelic frequency spectra that placed them together in clustering analysis as compared to the reference species populations. These results indicate that pines from the BS contact zone represent a highly hybridized population, and the presence of genetically pure individuals of P. sylvestris seems marginal, despite a number of macro-phenotypically trees deceptively resembling this species. Overall, pines from that area showed higher genetic similarity to P. uliginosa and P. uncinata, but not to P. mugo s.s. and P. sylvestris (Figure 2 and Figures S2 and S4). The predominance of hybrid individuals and their morphological differentiation could indicate either hybridization lasting for many generations with restricted inflow of genes of P. sylvestris, or more recent hybridization without competition of P. sylvestris individuals from the P. mugo complex. However, both hypotheses shall be verified in the special study applying more genetic markers at nuclear loci to reveal the proportion of genetic background from each parental species at much finer resolution.
Some discordance of needle characteristics and cpDNA markers in pines has been reported earlier for P. uncinata. Four individuals sampled as P. sylvestris appeared to be cryptic hybrids with cpDNA typical for the taxa of P. mugo and a set of needle characteristics that placed them at a marginal position of P. sylvestris individuals [31]. In the present biometric investigations, the majority of individuals determined during sampling as P. sylvestris did not resemble this species from the reference populations, neither in terms of needle (Figure 6A), nor in terms of cone characteristics (Figure 6B), falling between P. sylvestris and P. uliginosa reference populations (see also Figure S5). The only three individuals carrying P. sylvestris cpDNA also appeared close to P. sylvestris in the biometrical analyses on the needle and cone characteristics. During sampling, only one of them was classified as P. sylvestris based on the main phenotype traits used in the field of taxonomic classification of trees [15,16]. However, the remaining individuals classified as this species were intermediate between P. sylvestris and P. uliginosa, indicating a lack of correlation between habit form vs. needle and cone traits. Nevertheless, this observation has not been verified using statistical methods. The lack of correlation between morphological characteristics of hybrids has been reported earlier by Bobowicz et al. [89], but on restricted and relatively young material.
Interestingly, the individuals bearing cpDNA of P. sylvestris were close to the reference individuals of the species in terms of PCA and discrimination scatter plots for needle and cone characteristics (Figure 6A,B and Figure S6) showing that the results of the biometrical analyses are similar to those of cpDNA used for distinguishing between P. sylvestris and the P. mugo complex. In spite of that, the 64 individuals representing the BS mixed population formed one group situated between P. sylvestris and the taxa of the P. mugo complex, without a clear division between P. sylvestris and P. uliginosa, as detected by Jasińska et al. [31] for P. sylvestris and P. uncinata.
The presence of P. uliginosa × P. sylvestris hybrids bearing the phenotype of P. sylvestris and the lack of a vice versa morphological combination indicates a potential unidirectional gene flow from P. uliginosa to P. sylvestris, without phenotype expression of the former taxon (e.g., [30,31] and the literature cited herein). The putative gene flow from P. uliginosa (also from other taxa of the P. mugo complex) to P. sylvestris could be a specific evolutionary attribute, which allow the taxa of the P. mugo complex to diverge from P. sylvestris. The known invasions due to gene flow by pollen are a frequent phenomenon in nature. The gene inflow of Quercus petraea (Matt.) Liebl. into Q. robur L. by pollen is a good example of such a process [90,91,92], which also takes place in the populations of other wind-pollinated central-European trees ([93,94] and the literature cited herein). Examples of invasion by hybridization and gene flow were described also for other plant species [95,96,97], but only certain types of hybrids have adaptive advantages over other hybrid types and parental species in diverse environments. This hypothesis, however, shall be verified in the study on the genetic basis of hybrid adaptation, which requires in-depth investigation of polymorphisms at various genomic regions [98].
Predominance of individuals bearing cpDNA barcodes of the P. mugo complex suggests long-lasting hybridization in isolation from other populations of P. sylvestris, as a result of vegetation changes during the Holocene [44,45,48]. The occurrence of pines in the Sudetes was strongly reduced during the Boreal period of the Holocene. This was mainly due to the expansion of Picea, which pushed pines to specific sites, such as rocky ridges or rock tops and, consequently, involved the spatial isolation between their populations. Thus, in Błędne Skały we observe the result of hybridization between P. sylvestris and P. uliginosa without or with highly restricted gene inflow from other populations. Isolation together with unidirectional gene flow from the P. mugo complex to P. sylvestris, could explain the considerable presence of hybrid individuals there and the restricted number of specimens with a genome typical for the latter species.
Despite its hybrid origin, we detected a relatively low level of genetic variation in the mixed BS population in terms of mitochondrial and plastid DNA markers, lower than in the reference populations of P. mugo s.s., P. sylvestris and P. uliginosa [63,64,65,68,99,100]. Some discrepancy between the level of genetic variation at the genomes of different inheritance modes in pines suggests that a limited number of paternal trees contributed to the reproductive success of the population. However, this observation and scale of possible inbreeding would need some additional investigations including genetic variation of the progeny derived from seeds sampled in the population.

5. Conclusions

The investigated population consists mostly of hybrids, and the presence of pure P. sylvestris individuals seems marginal. Hybridization between P. sylvestris and the taxa of the P. mugo complex is predominantly unidirectional from the latter species to P. sylvestris. The macro-morphological characteristics of the latter parent tree are frequently conserved in hybrids, but we did not find correlation of those traits with morphological and anatomical characteristics of the needles and cone traits in our study.

Supplementary Materials

The following are available online at https://www.mdpi.com/1999-4907/11/10/1086/s1, Figure S1: Median Joining network of the mitotypes of pines from the mixed population in Błędne Skały and the reference four pine species populations based on 14 mtDNA regions, Figure S2: UPGMA tree based on Nei genetic distance at 14 mtDNA regions showing genetic relationships between defined groups of samples from the mixed pine population in Błędne Skały and reference pure species populations, Figure S3: ΔK and L(K) indices showing the most likely partitioning of the genetic variation at nSSR loci in Structure analysis, Figure S4: STRUCTURE results indicating three genetic groups of samples based on genetic variation at cpSSR loci (above) with optimal number of clusters detected using the method of Evanno’s ΔK among compared pine populations (below), Figure S5: Influence of needle (A) and cone (B) characteristics on the dispersion of samples of four sample groups from Błędne Skały and reference pine taxa, Figure S6: Results of discrimination analysis for mixed populations from Błędne Skały (BS), P. sylvestris (PS), P. mugo (PM), P. uliginosa (PUL) and P. uncinata (PUN) based on the characteristic of needles (A) and cones (B), Table S1: mtDNA haplotypes and their frequencies in the four groups of pines from Błędne Skały, Table S2: Frequency of the 20 most frequent cpDNA haplotypes (present in > 3 individuals) out of 302 detected ones in the four groups of pines from Błędne Skały, Table S3: Genetic variation at individual populations based on cpSSR loci, Table S4: Genetic variation of samples from the mixed pine population in Błędne Skały and reference populations of the four pine species based on nDNA SSR loci, Table S5: Basic summary statistics of analysed needle traits of pines from the Błędne Skały (BS) and of the taxa compared—P. sylvestris (S), P. mugo (M), P. uliginosa (UL) and P. uncinata (UN), Table S6: Basic summary statistics of analysed cone characteristics of pines from Błędne Skały (BS) and of the taxa compared—P. sylvestris (S), P. mugo (M), P. uliginosa (UL) and P. uncinata (UN), Table S7: Needle characteristics differing at p ≤ 0.01 and p ≤ 0.05 in the Tukey’s HSD and Kruskal–Wallis tests between subpopulations from Błędne Skały and compared populations of P. sylvestris, P. mugo, P. uliginosa and P. uncinata, Table S8: Cone characteristics differing at p ≤ 0.01 and p ≤ 0.05 in Tukey’s HSD and Kruskal–Wallis tests between subpopulations from Błędne Skały and compared populations of P. sylvestris, P. mugo, P. uliginosa and P. uncinata.

Author Contributions

Conceptualization, K.S. (Karolina Sobierajska), W.W and A.B.; methodology, W.W. and A.B.; formal molecular analysis, B.W., J.Z. and B.Ł.; formal biometric analysis, M.S., A.K.J., K.S. (Karolina Sobierajska), K.S. (Katarzyna Sękiewicz), K.B., K.M.; writing—original draft preparation, K.S. (Karolina Sobierajska), W.W., B.Ł. and A.B.; writing—review and editing, B.Ł. and W.W.; supervision, funding acquisition, W.W. and A.B. All authors have read and agreed to the published version of the manuscript.

Funding

This research was supported by Polish National Science Centre (UMO-2015/19/B/NZ9/00024) and the Stołowe Mountains National Park (ZBN/03/2020). We thank Director and management of the park for their help in sampling of plant material.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Rieseberg, L.H.; Ellstrand, N.C.; Arnold, M. What can molecular and morphological markers tell us about plant hybridization? Crit. Rev. Plant Sci. 1993, 12, 213–241. [Google Scholar] [CrossRef]
  2. Arnold, M.L. Natural Hybridization and Evolution; Oxford University Press: Oxford, UK, 1997. [Google Scholar]
  3. Arnold, M.L. Evolution through Genetic Exchange; Oxford University Press: Oxford, UK, 2007; pp. 1–252. [Google Scholar] [CrossRef]
  4. Mallet, J. Hybridization, ecological races and the nature of species: Empirical evidence for the ease of speciation. Philos. Trans. R Soc. B Biol. Sci. 2008, 363, 2971–2986. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  5. Soltis, P.S.; Soltis, D.E. The role of hybridization in plant speciation. Annu. Rev. Plant Biol. 2009, 60, 561–588. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  6. Abbott, R.; Albach, D.; Ansell, S.; Arntzen, J.W.; Baird, S.J.E.; Bierne, N.; Boughman, J.; Brelsford, A.; Buerkle, C.A.; Buggs, R.; et al. Hybridization and speciation. J. Evol. Biol. 2013, 26, 229–246. [Google Scholar] [CrossRef] [Green Version]
  7. Boratyński, A.; Boratyńska, K.; Lewandowski, A.; Gołąb, Z.; Kiciński, P. Evidence of the possibility of natural reciprocal crosses between Pinus sylvestris and P. uliginosa based on the phenology of reproductive organs. Flora 2003, 198, 377–388. [Google Scholar] [CrossRef]
  8. Stacy, E.A.; Paritosh, B.; Johnson, M.A.; Price, D.K. Incipient ecological speciation between successional varieties of a dominant tree involves intrinsic postzygotic isolating barriers. Ecol. Evol. 2017, 7, 2501–2512. [Google Scholar] [CrossRef]
  9. Gernandt, D.S.; López, G.G.; García, S.O.; Liston, A. Phylogeny and classification of Pinus. Taxon 2005, 54, 29–42. [Google Scholar] [CrossRef] [Green Version]
  10. Kormuťák, A.; Ostrolucká, M.; Vooková, B.; Preťová, A.; Fečková, M. Artificial hybridization of Pinus sylvestris L. and Pinus mugo Turra. Acta. Biol. Crac. Ser. Bot. 2005, 47, 129–134. [Google Scholar]
  11. Delgado, P.; Salas-Lizana, R.; Vazquez-Lobo, A.; Wegier, A.; Anzidei, M.; Alvarez-Buylla, E.R.; Vendramin, G.G.; Pinero, D. Introgressive hybridization in Pinus montezumae Lamb and Pinus pseudostrobus Lindl. (Pinaceae): Morphological and molecular (cpSSR) evidence. Int. J. Plant Sci. 2007, 168, 861–875. [Google Scholar] [CrossRef]
  12. Krajmerová, D.; Paule, L.; Zhelev, P.; Voleková, M.; Evtimov, I.; Gagov, V.; Gömöry, D. Natural hybridization in eastern-Mediterranean firs: The case of Abies borisii-regis. Plant Biosyst. 2016, 150, 1189–1199. [Google Scholar] [CrossRef]
  13. Wachowiak, W.; Palme, A.E.; Savolainen, O. Speciation history of three closely related pines Pinus mugo (T.), P. uliginosa (N.) and P. sylvestris (L.). Mol. Ecol. 2011, 20, 1729–1743. [Google Scholar] [CrossRef] [PubMed]
  14. Wachowiak, W.; Zaborowska, J.; Labiszak, B.; Perry, A.; Zucca, G.M.; Gonzalez-Martinez, S.C.; Cavers, S. Molecular signatures of divergence and selection in closely related pine taxa. Tree Genet. Genomes 2018, 14. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  15. Christensen, K.l. Taxonomic revision of the Pinus mugo complex and P. rhaetica (P. mugo sylvestris) (Pinaceae). Nord. J. Bot. 1987, 7, 383–408. [Google Scholar] [CrossRef]
  16. Businský, R.; Kirschner, J. Pinus mugo and P. uncinata as parents of hybrids a taxonomic and nomenclatural survey. Phyton 2010, 50, 27–57. [Google Scholar]
  17. Staszkiewicz, J.; Tyszkiewicz, M. Variability of the natural hybrids of Pinus sylvestris L. x Pinus mugo Turra (= P. x rotundata Link) in south-western Poland and in some selected localities of Bohemia and Moravia. Fragm. Flor. Geobot. 1972, 18, 173–191. [Google Scholar]
  18. Bobowicz, M.A. Differentiation of Pinus sylvestris L. and Pinus mugo Turra, pines from Bór na Czerwonem and from Zieleniec in traits of one- and two-year-old cones. Bull. Soc. Amis. Sci. Lett. Pozn. Ser. D Sci. Biol. 1988, 26, 99–108. [Google Scholar]
  19. Neet-Sarqueda, C. Genetic differentiation of Pinus sylvestris L. and Pinus mugo aggr. populations in Switzerland. 1994, 43, 207–214. [Google Scholar]
  20. Christensen, K.I.; Dar, G.H. A morphometric analysis of spontaneous and artificial hybrids of Pinus mugo ×sylvestris (Pinaceae). Nord. J. Bot. 1997, 17, 77–86. [Google Scholar] [CrossRef]
  21. Boratyńska, K.; Bobowicz, M.A. Pinus uncinata Ramond taxonomy based on needle characters. Plant. Syst. Evol. 2001, 227, 183–194. [Google Scholar] [CrossRef]
  22. Boratyńska, K.; Boratyński, A.; Lewandowski, A. Morphology of Pinus uliginosa (Pinaceae) needles from populations exposed to and isolated from the direct influence of Pinus sylvestris. Bot. J. Linn. Soc. 2003, 142, 83–91. [Google Scholar] [CrossRef] [Green Version]
  23. Boratyńska, K.; Boratyński, A. Taxonomic differences among closely related pines Pinus sylvestris, P. mugo, P. uncinata, P. rotundata and P. uliginosa as revealed in needle sclerenchyma cells. Flora 2007, 202, 555–569. [Google Scholar] [CrossRef]
  24. Boratyńska, K.; Jasińska, A.; Ciepłuch, E. Effect of tree age on needle morphology and anatomy of Pinus uliginosa and Pinus silvestris—Species-specific character separation during ontogenesis. Flora Morphol. Distrib. Funct. Ecol. Plants 2008, 203, 617–626. [Google Scholar] [CrossRef]
  25. Venturas, M.; García Álvarez, S.; Alcantra, M.F.; Collada, C.; Gil, L. Species selection for reforestations: What happens with historic local extinctions and habitat protection zones? A case study in the Cantabrian Range. Eur. J. Res. 2013, 132, 107–120. [Google Scholar] [CrossRef]
  26. Lewandowski, A.; Smoćko, J.; Boratynska, K.; Boratynski, A. Genetic differences between two Polish populations of Pinus uliginosa, compared to P. sylvestris and P. mugo. Dendrobiology 2002, 48, 51–57. [Google Scholar]
  27. Wachowiak, W.; Celinski, K.; Prus-Glowacki, W. Evidence of natural reciprocal hybridisation between Pinus uliginosa and P. sylvestris in the sympatric population of the species. Flora 2005, 200, 563–568. [Google Scholar] [CrossRef]
  28. Wachowiak, W.; Lewandowski, A.; Prus-Głowacki, W. Reciprocal controlled crosses between Pinus sylvestris and P. mugo verified by a species-specific cpDNA marker. J. Appl. Genet. 2005, 46, 41–43. [Google Scholar]
  29. Kormutak, A.; Vookova, B.; Manka, P.; Salaj, J.; Camek, V.; Gömöry, D. Abortive embryogenesis in hybrid swarm populations of Pinus sylvestris L. and Pinus mugo Turra. Trees 2008, 22, 657–662. [Google Scholar] [CrossRef]
  30. Wachowiak, W.; Prus-Glowacki, W. Hybridisation processes in sympatric populations of pines Pinus sylvestris L., P.mugo Turra and P.uliginosa Neumann. Plant Syst. Evol. 2008, 271, 29–40. [Google Scholar] [CrossRef]
  31. Jasińska, A.K.; Wachowiak, W.; Muchewicz, E.; Boratyńska, K.; Montserrat, J.M.; Boratyński, A. Cryptic hybrids between Pinus uncinata and P. sylvestris. Bot. J. Linn. Soc. 2010, 163, 473–485. [Google Scholar] [CrossRef]
  32. Wachowiak, W.; Cavers, S.; Żukowska, W.B. Interspecific gene flow and ecological selection in a pine (Pinus sp.) contact zone. Plant Syst. Evol. 2015, 301, 1643–1652. [Google Scholar] [CrossRef] [Green Version]
  33. Wachowiak, W.; Zukowska, W.B.; Wojkiewicz, B.; Cavers, S.; Litkowiec, M. Hybridization in contact zone between temperate European pine species. Tree Genet. Genomes 2016, 12. [Google Scholar] [CrossRef] [Green Version]
  34. Lewandowski, A.; Dering, M. Crossability between Pinus uliginosa and its putative parental species Pinus sylvestris and Pinus mugo. Silvae Genet. 2006, 55, 52–54. [Google Scholar] [CrossRef] [Green Version]
  35. Kormutak, A.; Galgoci, M.; Manka, P.; Koubova, M.; Jopcik, M.; Sukenikova, D.; Bolecek, P.; Gőmőry, D. Field-based artificial crossings indicate partial compatibility of reciprocal crosses between Pinus sylvestris and Pinus mugo and unexpected chloroplast DNA inheritance. Tree Genet. Genomes 2017, 13, 68. [Google Scholar] [CrossRef]
  36. Liston, A.; Parker-Defeniks, M.; Syring, J.; Willyard, A.; Cronn, R. Interspecific phylogenetic analysis enhances intraspecific phylogeographical inference: A case study in Pinus lambertiana. Mol. Ecol. 2007, 16, 3926–3937. [Google Scholar] [CrossRef]
  37. Senjo, M.; Kimura, K.; Watano, Y.; Ueda, K.; Shimizu, T. Extensive mitochondrial introgression from Pinus pumila to P. parviflora var. pentaphylla (Pinaceae). J. Plant Res. 1999, 112, 97–105. [Google Scholar] [CrossRef]
  38. Ma, X.-F.; Szmidt, A.E.; Wang, X.-R. Genetic structure and evolutionary history of a diploid hybrid pine pinus densata inferred from the nucleotide variation at seven gene loci. Mol. Biol. Evol. 2006, 23, 807–816. [Google Scholar] [CrossRef] [Green Version]
  39. Staffa, M. Słownik Geografii Turystycznej Sudetów: Góry Stołowe [Lexicon of the Turistic Geography of the Sudetes. 13. Stołowe Mountains]; PTTK Kraj: Kraków, Poland, 1992. [Google Scholar]
  40. Latocha, A.; Migoń, P. Geneza i przemiany krajobrazu kulturowego Gór Stołowych [Genezis and transformations of the cultural lanscape of the Stołowe Mountains]. In Góry Stołowe—Przyroda i Ludzie [The Stołowe Mountains—Nature and People]; Kabała, C., Ed.; Park Narodowy Gór Stołowych: Kudowa Zdrój, Poland, 2018; pp. 47–61. [Google Scholar]
  41. Latałowa, M.; Tobolski, K.; Nalepka, D. Pinus L. subgenus Pinus (subgen. Diploxylon (Koehne) Pilger). In Late Glacial and Holocene History of Vegetation in Poland Based on Isopollen Maps; Ralska-Jasiewiczowa, M., Ed.; W. Szafer Institute of Botany: Kraków; Poland, 2004; pp. 165–177. [Google Scholar]
  42. Treml, V.; Jankovská, V.; Petr, L. Holocene dynamics of the alpine timberline in the High Sudetes. Biologia 2008, 63, 73–80. [Google Scholar] [CrossRef]
  43. Malkiewicz, M.; Waroszewski, J.; Bojko, O.; Egli, M.; Kabala, C. Holocene vegetation history and soil development reflected in the lake sediments of the Karkonosze Mountains (Poland). Holocene 2016, 26, 890–905. [Google Scholar] [CrossRef] [Green Version]
  44. Marek, S. Rozwój Wielkiego Torfowiska Batorowskiego w świetle badań biostratygraficznych [Development of Wielkie Torfowisko Barowskie mire in the light of biostratygraphic investigations]. Szczeliniec 1988, 2, 49–88. [Google Scholar]
  45. Madeyska, E. The history of the Zieleniec Mire and the surrounding areas based on the palynological research. Monogr. Bot. 2005, 94, 145–157. [Google Scholar]
  46. Baranowska-Kącka, A. Historia roślinności Gór Stołowych w świetle badań palinologicznych Sudetów [The history of vegetation of the Stołowe Mts in the light of palynological studies]. In Przyroda Parku Narodowego Gór Stołowych [Nature of National Park of the Stołowe Mts]; Witkowski, A., Pokryszko, B.M., Ciężkowski, W., Eds.; National Park of the Stołowe Mountains: Kudowa Zdrój, Poland, 2013; pp. 130–136. [Google Scholar]
  47. Boratyński, A. Pinus uliginosa Neumann in the Błedne Skały reserve in Stolowe Mountains. Dendrobiology 1978, 23, 261–267. [Google Scholar]
  48. Glina, B.; Malkiewicz, M.; Mendyk, Ł.; Bogacz, A.; Woźniczka, P. Human-affected disturbances in vegetation cover and peatland development in the late Holocene recorded in shallow mountain peatlands (Central Sudetes, SW Poland). Boreas 2017, 46, 294–307. [Google Scholar] [CrossRef]
  49. Sarvas, R. Investigations on the flowering and seed crop of Pinus silvestris. Commun. Inst. Fenn. 1962, 53, 198. [Google Scholar]
  50. Sarvas, R. Investigations on the annual cycle of development of forest trees. Active period. Commun. Inst. Fenn. 1972, 76, 110. [Google Scholar]
  51. Poska, A.; Pidek, I.A. Pollen dispersal and deposition characteristics of Abies alba, Fagus sylvatica and Pinus sylvestris, Roztocze region (SE Poland). Veg. Hist. Archaeobot. 2010, 19, 91–101. [Google Scholar] [CrossRef]
  52. Robledo-Arnuncio, J.J.; Gil, L. Patterns of pollen dispersal in a small population of Pinus sylvestris L. revealed by total-exclusion paternity analysis. Heredity 2005, 94, 13–22. [Google Scholar] [CrossRef] [PubMed]
  53. Burczyk, J.; Chalupka, W. Flowering and cone production variability and its effect on parental balance in a Scots pine clonal seed orchard. Ann. Sci. 1997, 54, 129–144. [Google Scholar] [CrossRef] [Green Version]
  54. Marcysiak, K. Interpopulational variability of Pinus uncinata Ramond ex DC. in Lam. & DC (Pinaceae) cone characters. Dendrobiology 2004, 51, 43–51. [Google Scholar]
  55. Marcysiak, K.; Boratyńska, K.; Mazur, M. Variability of Pinus uliginosa cones from the peat-bog in Węgliniec. Dendrobiology 2003, 49, 43–47. [Google Scholar]
  56. Boratyńska, K.; Marcysiak, K.; Boratyński, A. Pinus mugo (Pinaceae) in the Abruzzi Mountains: High morphological variation in isolated populations. Bot. J. Linn. Soc. 2005, 147, 309–316. [Google Scholar] [CrossRef] [Green Version]
  57. Marcysiak, K.; Boratyński, A. Contribution to the taxonomy of Pinus uncinata (Pinaceae) based on cone characters. Plant Syst. Evol. 2007, 264, 57–73. [Google Scholar] [CrossRef]
  58. Sobierajska, K.; Boratynska, K. Variability of needle characters of Pinus mugo Turra populations in the Karkonosze Mountains in Poland. Dendrobiology 2008, 59, 41–49. [Google Scholar]
  59. Boratynska, K.; Lewandowska, D. Differences among three populations of Pinus uliginosa and their relation to P. sylvestris as expressed by the needle characters. Dendrobiology 2009, 61, 37–46. [Google Scholar]
  60. Sobierajska, K.; Boratyńska, K.; Marcysiak, K. Variation of cone characters in Pinus mugo (Pinaceae) populations in the Giant Mountains (Karkonosze, Sudetes). Dendrobiology 2010, 63, 33–41. [Google Scholar]
  61. Jasińska, A.; Boratyńska, K.; Dering, M.; Sobierajska, K.; Ok, T.; Romo Diez, A.M.; Boratyński, A. Distance between south-European and south-west Asiatic refugial areas involved morphological differentiation: Pinus sylvestris case study. Plant Syst. Evol. 2014, 300, 1487–1502. [Google Scholar] [CrossRef] [Green Version]
  62. Dzialuk, A.; Boratyńska, K.; Romo, A.; Boratyński, A. Taxonomic and geographic variation of the Pinus mugo complex on chloroplast microsatellite markers. Syst. Biodivers. 2016, 15, 464–479. [Google Scholar] [CrossRef] [Green Version]
  63. Wachowiak, W.; Boratynska, K.; Cavers, S. Geographical patterns of nucleotide diversity and population differentiation in three closely related European pine species in the Pinus mugo complex. Bot. J. Linn. Soc. 2013, 172, 225–238. [Google Scholar] [CrossRef] [Green Version]
  64. Łabiszak, B.; Zaborowska, J.; Wachowiak, W. Patterns of mtDNA variation reveal complex evolutionary history of relict and endangered peat bog pine (Pinus uliginosa). AoB Plants 2019, 11. [Google Scholar] [CrossRef]
  65. Boratyńska, K.; Dzialuk, A.; Lewandowski, A.; Marcysiak, K.; Jasińska, A.; Sobierajska, K.; Tomaszewski, D.; Burczyk, J.; Boratyński, A. Geographic distribution of quantitative traits variation and genetic variability in natural populations of Pinus mugo in Central Europe. Dendrobiology 2014, 72, 65–84. [Google Scholar] [CrossRef] [Green Version]
  66. Donnelly, K.; Cottrell, J.; Ennos, R.A.; Vendramin, G.G.; A’Hara, S.; King, S.; Perry, A.; Wachowiak, W.; Cavers, S. Reconstructing the plant mitochondrial genome for marker discovery: A case study using Pinus. Mol. Ecol. Resour. 2017, 17, 943–954. [Google Scholar] [CrossRef] [Green Version]
  67. Soranzo, N.; Alia, R.; Provan, J.; Powell, W. Patterns of variation at a mitochondrial sequence-tagged-site locus provides new insights into the postglacial history of European Pinus sylvestris populations. Mol. Ecol. 2000, 9. [Google Scholar] [CrossRef] [PubMed]
  68. Zaborowska, J.; Łabiszak, B.; Wachowiak, W. Population history of European mountain pines Pinus mugo and Pinus uncinata revealed by mitochondrial DNA markers. J. Syst. Evol. 2019. [Google Scholar] [CrossRef]
  69. Żukowska, W.B.; Boratyńska, K.; Wachowiak, W. Comparison of range-wide chloroplast microsatellite and needle trait variation patterns in Pinus mugo Turra (dwarf mountain pine). iForest 2017, 10, e1–e9. [Google Scholar] [CrossRef] [Green Version]
  70. Wachowiak, W.; Lesniewicz, K.; Odrzykoski, I.; Augustyniak, H.; Prus-Głowacki, W. Species specific cpDNA markers useful for studies on the hybridisation between Pinus mugo and P. sylvestris. Acta Soc. Bot. Pol. 2000, 69, 273–276. [Google Scholar] [CrossRef] [Green Version]
  71. Żukowska, W.B.; Wachowiak, W. Nuclear microsatellite markers reveal the low genetic structure of Pinus mugo Turra (dwarf mountain pine) populations in Europe. Plant Syst. Evol. 2017, 303, 641–651. [Google Scholar] [CrossRef]
  72. Bandelt, H.-J.; Forster, P.; Röhl, A. Median-joining networks for inferring intraspecific phylogenies. Mol. Biol. Evol. 1999, 16, 37–48. [Google Scholar] [CrossRef]
  73. Leigh, J.W.; Bryant, D. Popart: Full-feature software for haplotype network construction. Methods Ecol. Evol. 2015, 6, 1110–1116. [Google Scholar] [CrossRef]
  74. Librado, P.; Rozas, J. DnaSP v5: A software for comprehensive analysis of DNA polymorphism data. Bioinformatics 2009, 25, 1451–1452. [Google Scholar] [CrossRef] [Green Version]
  75. Kumar, S.; Stecher, G.; Tamura, K. MEGA7: Molecular evolutionary genetics analysis version 7.0 for bigger datasets. Mol. Biol. Evol. 2016, 33, 1870–1874. [Google Scholar] [CrossRef] [Green Version]
  76. Excoffier, L.; Lischer, H.E.L. Arlequin suite ver 3.5: A new series of programs to perform population genetics analyses under Linux and Windows. Mol. Ecol. Resour. 2010, 10, 564–567. [Google Scholar] [CrossRef]
  77. Peakall, R.; Smouse, P.E. GenAlEx 6.5: Genetic analysis in Excel. Population genetic software for teaching and research-an update. Bioinformatics 2012, 28, 2537–2539. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  78. Pritchard, J.K.; Stephens, M.; Donnelly, P. Inference of population structure using multilocus genotype data. Genetics 2000, 155, 945–959. [Google Scholar] [PubMed]
  79. Falush, D.; Stephens, M.; Pritchard, J.K. Inference of population structure using multilocus genotype data: Dominant markers and null alleles. Mol. Ecol. Notes 2007, 7, 574–578. [Google Scholar] [CrossRef] [PubMed]
  80. Hubisz, M.J.; Falush, D.; Stephens, M.; Pritchard, J.K. Inferring weak population structure with the assistance of sample group information. Mol. Ecol. Resour. 2009, 9, 1322–1332. [Google Scholar] [CrossRef] [Green Version]
  81. Evanno, G.; Regnaut, S.; Goudet, J. Detecting the number of clusters of individuals using the software STRUCTURE: A simulation study. Mol. Ecol. 2005, 14, 2611–2620. [Google Scholar] [CrossRef] [Green Version]
  82. Earl, D.A.; Vonholdt, B.M. STRUCTURE HARVESTER: A website and program for visualizing STRUCTURE output and implementing the Evanno method. Conserv. Genet. Resour. 2012, 4, 359–361. [Google Scholar] [CrossRef]
  83. Boratyńska, K.; Tomaszewski, D.; Montserrat Martí, J.; Marek, S.; Boratyński, A. Taxonomic position of Pinus ceciliae (Pinaceae) endemic for Balearic Islands as revealed on needle characteristics. Dendrobiology 2019, 82, 8–16. [Google Scholar] [CrossRef]
  84. Staszkiewicz, J. Morphological variation of leaves, cones and seeds. In Scots Pine Biology; Białobok, S., Boratyński, A., Bugała, W., Eds.; Instytut Dendrologii Polskiej Akademii Nauk: Poznań—Kórnik, Poland, 1993; pp. 33–44. [Google Scholar]
  85. Zar, J.H. Biostatistical Analysis, 5th ed.; Prentice-Hall/Pearson: Upper Saddle River, NJ, USA, 2007. [Google Scholar]
  86. Sokal, R.R.; Rohlf, F.J. Biometry, 4th ed.; W.H. FREEMAN & CO LTD: New York, NY, USA, 2012. [Google Scholar]
  87. Watala, C. Biostatystyka: Wykorzystanie Metod Statystycznych w Pracy Badawczej w Naukach Biomedycznych [Biostatistics: Using Statistic Methods in Biomedical Scientific Investigations]; Alfamedica Press: Bielsko Biała, Poland, 2012. [Google Scholar]
  88. Morrison, D.F. Multivariate Statistical Methods; Thomson/Brooks/Cole: Belmont, CA, USA, 2005. [Google Scholar]
  89. Bobowicz, M.; Danielewicz, W. Isoenzymatic variability in progeny of Pinus mugo Turra x Pinus sylvestris L. hybrids from Bór na Czerwonem, in experimental culture. Acta Soc. Bot. Pol. 2000, 69, 137–144. [Google Scholar] [CrossRef]
  90. Kremer, A.; Dupouey, J.-L.; Deans, J.; Cottrell, J.; Csaikl, U.; Finkeldey, R.; Espinel, S.; Jensen, J.; Kleinschmit, J.; Dam, B.C.; et al. Leaf morphological differentiation between Quercus robur and Quercus petraea is stable across western European mixed oak stands. Ann. For. Sci. 2002, 59, 777–787. [Google Scholar] [CrossRef]
  91. Petit, R.J.; Bodénès, C.; Ducousso, A.; Roussel, G.; Kremer, A. Hybridization as a mechanism of invasion in oaks. New Phytol. 2004, 161, 151–164. [Google Scholar] [CrossRef] [Green Version]
  92. Boratyński, A.; Marcysiak, K.; Jasińska, A.K.; Iszkuło, G. Interrelations among con-generic and co-occurring tree species: Asymmetric hybridization and the high success of Quercus petraea (Matt.) Liebl. regeneration in mixed Q. petraea Q. robur L. stands. Pol. J. Ecol. 2010, 58, 273–283. [Google Scholar]
  93. Petit, R.J. Biological invasions at the gene level. Divers. Distrib. 2004, 10, 159–165. [Google Scholar] [CrossRef]
  94. Petit, R.; Bialozyt, R.; Garnier-Géré, P.; Hampe, A. Ecology and genetics of tree invasions: From recent introductions to Quaternary migrations. Ecol. Manag. 2004, 197, 117–137. [Google Scholar] [CrossRef]
  95. Martinsen, G.D.; Whitham, T.G.; Turek, R.J.; Keim, P. Hybrid populations selectively filter gene introgression between species. Evolution 2001, 55, 1325–1335. [Google Scholar] [CrossRef]
  96. Fay, M.F.; Cowan, R.S.; Simpson, D.A. Hybridisation between Schoenoplectus tabernaemontani and S. triqueter (Cyperaceae) in the British Isles. Watsonia 2003, 24, 433–442. [Google Scholar]
  97. Álvarez, I.; Wendel, J.F. Cryptic Interspecific Introgression and Genetic Differentiation within Gossypium Aridum (Malvaceae) and Its Relatives. Evolution 2006, 60, 505–517. [Google Scholar] [CrossRef] [Green Version]
  98. Wachowiak, W.; Trivedi, U.; Perry, A.; Cavers, S. Comparative transcriptomics of a complex of four European pine species. BMC Genom. 2015, 16, 234. [Google Scholar] [CrossRef] [Green Version]
  99. Dzialuk, A.; Muchewicz, E.; Boratyński, A.; Montserrat, J.M.; Boratyńska, K.; Burczyk, J. Genetic variation of Pinus uncinata (Pinaceae) in the Pyrenees determined with cpSSR markers. Plant Syst. Evol. 2009, 277, 197–205. [Google Scholar] [CrossRef]
  100. Dzialuk, A.; Boratyński, A.; Boratyńska, K.; Burczyk, J. Geographic patterns of genetic diversity of Pinus mugo (Pinaceae) in Central European mountains. Dendrobiology 2012, 68, 31–41. [Google Scholar]
Forests 11 01086 g001 550
Figure 1. Geographic position of Błędne Skały mixed population (BS) and samples of the reference species—Pinus sylvestris (S), P. mugo s.s. (M), P. uliginosa (UL) and P. uncinata (UN); population acronyms as in Table 1.
Figure 1. Geographic position of Błędne Skały mixed population (BS) and samples of the reference species—Pinus sylvestris (S), P. mugo s.s. (M), P. uliginosa (UL) and P. uncinata (UN); population acronyms as in Table 1.
Forests 11 01086 g001
Forests 11 01086 g002 550
Figure 2. Median Joining network of the mitotypes of pines from the mixed population in Błędne Skały (BS) and the reference four pine species populations based on 14 mtDNA regions. Circle sizes represent relative frequencies of haplotypes in a combined sample, lines between two circles represent singular mutation events, and black dots substitute absent haplotypes (population acronyms as in Table 1).
Figure 2. Median Joining network of the mitotypes of pines from the mixed population in Błędne Skały (BS) and the reference four pine species populations based on 14 mtDNA regions. Circle sizes represent relative frequencies of haplotypes in a combined sample, lines between two circles represent singular mutation events, and black dots substitute absent haplotypes (population acronyms as in Table 1).
Forests 11 01086 g002
Forests 11 01086 g003 550
Figure 3. Principal coordinates analysis (PCoA) based on Nei genetic distance at nSSR loci showing genetic relationships between groups of samples defined at the mixed pine population from Błędne Skały (BS) and reference populations of the four pine species (population acronyms as in Table 1).
Figure 3. Principal coordinates analysis (PCoA) based on Nei genetic distance at nSSR loci showing genetic relationships between groups of samples defined at the mixed pine population from Błędne Skały (BS) and reference populations of the four pine species (population acronyms as in Table 1).
Forests 11 01086 g003
Forests 11 01086 g004 550
Figure 4. STRUCTURE results indicating five groups of samples based on genetic variation at nSSR loci of the mixed pine population from Błędne Skały (BS) and reference populations of the four pine species (population acronyms as in Table 1).
Figure 4. STRUCTURE results indicating five groups of samples based on genetic variation at nSSR loci of the mixed pine population from Błędne Skały (BS) and reference populations of the four pine species (population acronyms as in Table 1).
Forests 11 01086 g004
Forests 11 01086 g005 550
Figure 5. Principal component analysis (PCA) plot showing relationships between groups of samples defined at the mixed pine population in Błędne Skały as P. sylvestris (BS_S), P. uliginosa (BS_M), hybrids (BS_H), undetermined (BS_N) and reference populations of the four pine species (population acronyms as in Table 1): (A) based on needle traits; (B) based on cone traits.
Figure 5. Principal component analysis (PCA) plot showing relationships between groups of samples defined at the mixed pine population in Błędne Skały as P. sylvestris (BS_S), P. uliginosa (BS_M), hybrids (BS_H), undetermined (BS_N) and reference populations of the four pine species (population acronyms as in Table 1): (A) based on needle traits; (B) based on cone traits.
Forests 11 01086 g005
Forests 11 01086 g006 550
Figure 6. Principal component analysis PCA plot showing relationships between individuals from four groups defined at the mixed pine population in Błędne Skały (marked in circles) and reference species, including P. sylvestris (S), P. mugo (M), P. uliginosa (UL) and P. uncinata (UN); (numbers 17, 70 and 76 indicate individuals with cpDNA haplotype characteristic for P. sylvestris): (A) based on needle traits; (B) based on cone traits.
Figure 6. Principal component analysis PCA plot showing relationships between individuals from four groups defined at the mixed pine population in Błędne Skały (marked in circles) and reference species, including P. sylvestris (S), P. mugo (M), P. uliginosa (UL) and P. uncinata (UN); (numbers 17, 70 and 76 indicate individuals with cpDNA haplotype characteristic for P. sylvestris): (A) based on needle traits; (B) based on cone traits.
Forests 11 01086 g006
Table
Table 1. Plant material used in the study.
Table 1. Plant material used in the study.
SpeciesCodeLocationLatitude
N		Longitude
E		Altitude
(m)		Nb (needles/cones)Ng
Mixed
BS		BSBłedne Skały, Stołowe Mountains (Poland)50.4816.2985064 (320/315)64
Pinus sylvestris
S		S1St. Miguel d’Engolasters (Andorra)42.511.561550-32
S2Divčibare Mts. (Serbia)44.1019.99957-26
S3Tatra National Park (Poland)49.2719.831000-12
S4Beskid (Poland)49.4020.82600-45
S5Pieniny National Park (Poland)49.4220.36800-21
S6Stołowe Mountains (Poland)50.4316.2390030 (300/0)
S7Sudetes, Karkonosze (Poland)50.8315.64600-20
S8Silesian Lowland, Węgliniec (Poland)51.2815.2317034 (340/46)
S9Liski Nature Reserve (Poland)51.9322.82150-11
S10Bory Tucholskie Forest, (Poland)53.3217.949050 (0/50)
S11Sosny Taborskie Reserve (Poland)53.7820.02115-13
S12Joutsa (Finland)61.7626.10125-25
Pinus mugo
M		M1Dinaric Alps (Bosnia and Herzegovina)43.7218.242020-25
M2Carpathians, Fagaraş (Romania)45.6124.601900-25
M3Alps, Carnic Alps, Pramollo (Italy)46.5613.281600-25
M4Ammergau Alps (Germany)47.5310.921870-25
M5Tatra National Park (Poland)49.1220.05170030 (300/330)12
M6Tatra National Park, Grześ (Poland)49.2319.771620-25
M7Sudetes (Poland)50.7815.58135033 (330/310)30
M8Sudetes, Śnieżka (Poland)50.7415.741435-30
Pinus uliginosa
UL		UL1Alps, Mittelwalde (Germany)47.4811.27860-20
UL2Zieleniec Nature Reserve (Poland)49.1320.0370030 (300/50)
UL3Large Batorów Peatland (Poland)50.4516.3875050 (500/50)31
UL4Węgliniec Nature Reserve (Poland)51.2915.2319052 (520/230)30
Pinus uncinata
UN		UN1Vall de Nuria (Spain)42.2202.10220028 (280/50)
UN2Vall de Ransol (Andorra)42.3801.37200032 (320/49)
UN3Coll de Jau (France)42.692.251500-24
UN4Belagua (Spain)42.97-0.771760-24
UN5La Trapa near Jaca (Spain)42.69-0.541720-22
UN6Col de la Croix-Morand (France)45.602.851400-14
Nb—number of trees used in biometric analysis; Ng—number of trees used in genetic analysis.
Table 2. Genetic variation of samples from mixed pine population in Błędne Skały and four reference pine species based on: (A) chloroplast microsatellite (cpSSR) loci; (B) mitochondrial (mtDNA) markers; (C) nuclear microsatellite (nSSR) loci.
Table
(A)
(A)
Sample NameMeasures of Genetic Variation
NAhAeAPARHdD2sh
Błędne Skały641913157.560.88132.1
P. mugo105100994695.870.99913.5
P. sylvestris12811811846109.280.9995.6
P. uncinata4833323220.570.9727.2
P. uliginosa8140322822.700.9688.2
Table
(B)
(B)
Sample NameMeasures of Genetic Variation
NSSpdHnHpHd (SD)uH
Błędne Skały64901.83740.487 (0.071)0.495
P. mugo341024.03720.786 (0.034)0.810
P. sylvestris61903.311040.796 (0.038)0.809
P. uncinata44511.06760.687 (0.055)0.703
P. uliginosa1011111.5423160.904 (0.014)0.913
Table
(C)
(C)
Sample NameMeasures of Genetic Variation
NNaNeHoHeF
Błędne Skały634.782.340.4250.4670.046
P. mugo1054.782.060.3860.4600.142
P. sylvestris1287.332.370.4160.4510.076
P. uncinata483.892.170.3720.4100.051
P. uliginosa814.782.090.3930.4380.101
N—sample size, Ah—mean number of haplotypes, Ae—mean effective number of haplotypes, AP—number of private haplotypes, AR—haplotype richness, Hd—haplotype diversity, D2sh—mean genetic distance of individuals within populations, S—number of segregating sites, SP—number of private polymorphisms, d—average number of nucleotide differences between two sequences, Hn—number of haplotypes, Hp—number of private haplotypes, uH—unbiased gene diversity, Na—mean number of alleles per locus, Ne—mean effective number of alleles per locus, Ho—mean observed heterozygosity, He—mean expected heterozygosity, F—mean fixation index.

Share and Cite

MDPI and ACS Style

Sobierajska, K.; Wachowiak, W.; Zaborowska, J.; Łabiszak, B.; Wójkiewicz, B.; Sękiewicz, M.; Jasińska, A.K.; Sękiewicz, K.; Boratyńska, K.; Marcysiak, K.; et al. Genetic Consequences of Hybridization in Relict Isolated Trees Pinus sylvestris and the Pinus mugo Complex. Forests 2020, 11, 1086. https://doi.org/10.3390/f11101086

AMA Style

Sobierajska K, Wachowiak W, Zaborowska J, Łabiszak B, Wójkiewicz B, Sękiewicz M, Jasińska AK, Sękiewicz K, Boratyńska K, Marcysiak K, et al. Genetic Consequences of Hybridization in Relict Isolated Trees Pinus sylvestris and the Pinus mugo Complex. Forests. 2020; 11(10):1086. https://doi.org/10.3390/f11101086

Chicago/Turabian Style

Sobierajska, Karolina, Witold Wachowiak, Julia Zaborowska, Bartosz Łabiszak, Błażej Wójkiewicz, Maciej Sękiewicz, Anna K. Jasińska, Katarzyna Sękiewicz, Krystyna Boratyńska, Katarzyna Marcysiak, and et al. 2020. "Genetic Consequences of Hybridization in Relict Isolated Trees Pinus sylvestris and the Pinus mugo Complex" Forests 11, no. 10: 1086. https://doi.org/10.3390/f11101086

APA Style

Sobierajska, K., Wachowiak, W., Zaborowska, J., Łabiszak, B., Wójkiewicz, B., Sękiewicz, M., Jasińska, A. K., Sękiewicz, K., Boratyńska, K., Marcysiak, K., & Boratyński, A. (2020). Genetic Consequences of Hybridization in Relict Isolated Trees Pinus sylvestris and the Pinus mugo Complex. Forests, 11(10), 1086. https://doi.org/10.3390/f11101086

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop