The evolution of spring fen ecotypes in Rhinanthus: genetic evidence for parallel origins in Scandinavia after the last ice age

Locally adapted ecotypes can constitute an important part of the biodiversity, especially in young floras with few endemic species. However, the origins, distinctness and conservation value of many ecotypes remain uncertain because genetic data are lacking or no common-garden study has been carried out. In the present study, we evaluated the distinctness and genetic structure of a phenotypically deviating morph of Rhinanthus angustifolius, growing in calcareous spring fens on the Baltic island of Gotland. Using data from a common-garden experiment and analyses of nuclear microsatellite variation, we compared fen populations on Gotland with conspecific populations from habitats more typical of the study species. We also included the fen specialist R. osiliensis from the Baltic island of Saaremaa in the molecular analyses to make further inferences about the origin of the Gotlandic fen morph. Our data indicate that the Gotlandic fen populations constitute a phenotypically and genetically distinct ecotype that most likely has evolved at least two times on Gotland after the last ice age. In congruence with previous studies, we also infer that fen ecotypes have evolved independently on Gotland and Saaremaa. We propose a varietal status for the Gotlandic fen ecotype and give recommendations for the conservation of this taxon.


Introduction
During the last glacial period, which lasted until about 15.000 ybp, Scandinavia was covered by the thick Weichselian ice sheet (Björck 1995). As a consequence, the presentday flora is young with few endemic taxa at the species level (Jonsell 1988(Jonsell , 2004. Many plant species have, however, adapted to local ecological conditions through the formation of phenotypically distinct ecotypes (Turesson 1922).
Ecotypes constitute an important part of Scandinavian biodiversity and may be unique, not only to specific habitats but also to the region as a whole (Jonsell 1988). For conservation authorities, it is important to consider ecotypes and other types of habitat-correlated variation (e.g. clines; Gregor 1938) in order to conserve the genetic diversity and long-term adaptive potential of species that occur in Scandinavia (Lundqvist et al. 2007). This work, however, may be hampered by the lack of consensus regarding the taxonomic treatment of ecotypes, which can be variously treated as subspecies or varieties, or remain unnamed, depending on local traditions and other factors (Hamilton and Reichard 1992;Jonsell 2004). To assess the origins, distinctness and conservation value of ecotypes, it is appropriate to use a combination of molecular-genetic analyses (Lowe et al. 2004;Stuessy et al. 2014) and common-garden experiments (Heslop-Harrison 1964), and to apply commonly accepted principles to the taxonomic treatment of the entities recognized.
Species of the annual, hemiparasitic plant genus Rhinanthus (Orobanchaceae), especially the two most widely distributed ones, R. minor L. and R. angustifolius C.C.Gmel., are highly polymorphic and renowned for their capacity to undergo ecotypic differentiation in response to differences in selection pressure (ter Borg 1972;Zopfi 1993aZopfi , b, 1995. Several species have independently and recurrently evolved seasonal ecotypes, differing in flowering time and related features, but considerable variation has also been recorded in characters unrelated to phenology such as leaf shape, hairiness and flower morphology (de Soó 1929;de Soó and Webb 1972). Thus, patterns of variation in Rhinanthus species tend to be complex and reticulate, which has resulted in a lack of consensus regarding the number and taxonomic level of described taxa. Ecotypes of R. angustifolius, for example, have variously been recognized as varieties, subspecies or even species (von Sterneck 1901;de Soó 1929;Mizianty 1978). The taxonomy is complicated by extensive hybridization between species (Kwak 1980;Ducarme and Wesselingh 2005) and the existence of intermediate populations that inter-connect taxonomically recognized ecotypes (Karlsson 1974). Besides, several of the characters commonly used in the taxonomic literature are known to be plastic, especially in relation to the identity of the host species (Jonstrup et al. 2016). Hence, there is considerable potential for environmental factors, such as differences in host plant composition, to mask or obscure genetically determined differences. Molecular studies of R. minor, R. angustifolius and R. alectorolophus (Scop.) Polich attest to the complex patterns expressed in previous taxonomic treatments of species in this genus (Houston and Wolff 2012;Pleines et al. 2013;Jonstrup et al. 2018).
The species at focus in the present study, R. angustifolius, is normally found in open dry grasslands such as hay meadows, coastal meadows and road verges, but putative conspecifics also occur in calcareous spring fens on the Baltic island of Gotland (SE Sweden) (Marklund 1982;Martinsson 1997;Johansson et al. 2016). The fen plants are similar to other plants of R. angustifolius in floral characters and leaf morphology but have a later flowering time, produce more nodes and fewer and shorter branches and possess a distinctive red colouration due to anthocyanin pigmentation, unlike the green colour normally found in R. angustifolius (Lindell 2006). On the island of Saaremaa (W Estonia), just 200 km away from Gotland, yet another form of Rhinanthus is found in calcareous spring fens. This form has been named R. osiliensis (Ronniger & Saarsoo) Vassilcz. and has long been considered as a local endemic for Saaremaa by Estonian botanists (Ronniger 1934). Fen plants on Gotland share many features with R. osiliensis, but possess fewer and shorter glandular hairs than those of the Estonian taxon (Lindell 2006). At present, the few fen populations on Gotland are provisionally treated as R. osiliensis and listed as near threatened (NT) in the Swedish Redlist (ArtDatabanken 2015).
A recent molecular study has revealed that R. osiliensis and the Gotland fen populations constitute distinct evolutionary units with putatively independent origins (Talve et al. 2014a, b). However, the distinctness and conservation value of the Gotland fen plants remain uncertain, as does their genetic relationship to the more typical R. angustifolius phenotype. In the study presented here, we analysed nuclear microsatellite variation as well as phenotypic data from a common-garden experiment to answer the following questions: 1. Is there a phenotypically and genetically distinct entity of R. angustifolius related to calcareous spring fens on Gotland? 2. Are the fen populations on Gotland related to R. osiliensis on Saaremaa and how do they relate to the widespread R. angustifolius? 3. If there is a distinct entity on Gotland, should it be treated as a taxon of its own and does it constitute one or several conservation units at the genetic level?

Study species
Rhinanthus angustifolius is a summer annual hemiparasitic herb that grows in various types of open grassland in western Eurasia. Plants of R. angustifolius are 20-50 cm tall with a glabrous to moderately hairy stem, opposite simple leaves and a green to reddish colour. Bracts and calyces are glabrous, but the margins of the calyces may be scabrid or puberulent. The yellow, zygomorphic flowers are placed in a terminal raceme on the central stem but can also be produced on primary and secondary branches (ter Borg 1972;de Soó and Webb 1972). The flowering period spans from May to September, with great overlap between populations (Jonstrup et al. 2016). Individual flowering date is partly determined by the number of nodes below the first flower; as a result, early flowering individuals generally have the fewest nodes (Wesselingh 2016). The flowers are mainly pollinated by bumblebees and produce large, winged seeds in flat capsules (Kwak 1978;Natalis and Wesselingh 2012). The heavy seeds fall down close to the maternal plant when shaken by wind or other agents, but long-distance dispersal by animals or hay transfer is also possible (ter Borg 1972). The outcrossing rate within populations has been estimated to vary between 0.32 and 1.00 (mean 0.76; Ducarme and Wesselingh 2013). Plants of R. angustifolius have reduced root systems and attach to roots of other plants to support their requirement of water and essential nutrients (ter Borg 1972). The life cycle can be completed without attaching to a host, but then resulting in small plants with few flowers and seeds (A. Jonstrup, pers. obs.).

Seed collection, cultivation and phenotypic scoring
The Gotland fen phenotype has been reported from approximately six sites in two main areas on Gotland (Johansson et al. 2016). In August 2011, we collected seeds from three of these populations as well as from three geographically adjacent populations that conform to the more common and typical R. angustifolius phenotype (hereafter referred to as the fen morph and the common morph, respectively; Fig. 1, Table 1). In this study, we have chosen not to distinguish between two previously described taxa of R. angustifolius (ssp. angustifolius and ssp. grandiflorus (Wallr.) D.A. Webb) since results from common-garden experiments and molecular analyses contradict this division (Jonstrup et al. 2016(Jonstrup et al. , 2018. The two southern populations of the fen morph (SO2, SO3) inhabit fens grazed by Gotland ponies, while the northern one (SO1) inhabits a previously drained fen with no current management. Two of the common morph populations represent road verges (A11, A13) and the third an unmanaged coastal grassland (A12, Table 1). Capsules with ripe seeds were collected from 20 plants in each of the six study populations. To avoid collecting from closely related genotypes, we chose plants that grew > 2 m apart. We stored the seeds in separate paper bags at room temperature until they were sown in late November-early December the same year.
Seeds from each maternal plant were sown in three separate pots, filled with standard soil, consisting of peat (93%), clay (4%) and sand (3%), and the following contents of  Table 1 for full population names) essential macronutrients (nitrogen: 182 g/m 3 ; phosphorus: 91 g/m 3 and potassium: 194 g/m 3 ; pH 5.5-6.5). The pots were planted with one of three host species: Festuca rubra (Poaceae), Lolium perenne (Poaceae) and Trifolium repens (Fabaceae) (hereafter referred to as Festuca, Lolium and Trifolium, respectively). Seeds of the host plants were obtained from a commercial seed supplier (Weibulls Horto) and had been sown earlier during the autumn of 2011. Because many Lolium and Trifolium plants died during the winter, additional seeds of these species were sown in the spring of 2012. Pots were placed randomly within and between four outdoor cold frames (blocks) in a sunny part of an experimental garden at Lund University. The pots were watered, but no extra fertilizer was added.
The Rhinanthus seeds started to germinate in late March 2012. Excess seedlings were continuously removed until only two Rhinanthus plants remained in each pot. In the summer, we scored the Rhinanthus plants for phenotypic characters deemed taxonomically important according to the previous literature (Ronniger 1934;de Soó and Webb 1972;Oja and Talve 2012) (Table 2). For scoring of hairiness, a single flower with its adjoining bract was collected and stored in 70% ethanol until measurements were taken under a stereomicroscope. The density of glandular hairs on the calyx margin was determined by counting all glandular hairs within a single ocular field at 4× magnification, centred along the upper margin at the widest part of the calyx. The mean length of the three longest glandular hairs (expressed in units of gland widths) in the ocular field was also determined, as was the corresponding measure for the bract (Table 2). During the course of the experiment, we had problems with fungal infections and aphids, which resulted in some deformation of flowers and inflorescences. To account for these factors, we estimated the rate of deformation for each plant (scored visually on a scale of 0-2) to serve as a covariate in some analyses (see below).
Voucher specimens, preserved from plants cultivated in the common garden (see Table 1), have been deposited in the Botanical Museum, Lund University (LD).

Statistical analyses of phenotypic data
We performed a combination of univariate and multivariate analyses in SPSS version 23 (IBM) to evaluate the phenotypic distinctness of the fen morph and the common morph on Gotland. Because of the unbalanced experimental design resulting from missing data, it was necessary to pool data across maternal seed families from the same population. Based on visual inspections of residuals, the error variation for quantitative traits was approximately normally distributed; thus, no transformation was deemed necessary. The overall pattern of variation between populations and environments was assessed in a multivariate analysis of variance (MANOVA) with flowering start and measures of overall plant morphology as dependent variables; morph, population (nested within morph), host species, block and the morph × host species and population × host species interactions as independent variables; and rate of infection as covariate. Population was treated as fixed (as all the other factors) since a large part of the existing fen populations were included in the analysis, i.e. the studied populations could not be considered as a strict sample. Hair characters were analysed separately because of large gaps in the data set.
To visualize patterns of variation between morphs and populations, we performed a discriminant function analysis (DFA; also known as canonical variates analysis, Dunn and Everitt 1982), based on data pooled across blocks, with each population grown on a particular host species being a group and individual plants serving as within-group replicates. A DFA creates a reduced set of discriminant functions (DFs) that maximize the separation between pre-determined groups relative to the variation within groups (Dunn and Everitt 1982).
As a complementary approach, we assessed the phenotypic distinctness of the fen morph and the common morph by pooling quantitative data over populations, blocks and host environments, and subjecting the data to a DFA with morphs as groups. This approach mimics a natural situation in which the measured Rhinanthus plants represent a mix of sites and potentially co-occurring microhabitats (e.g. different host species) (Jonstrup et al. 2016). The distinctness of morphs was quantified as the percentage of plants classified into the original morph based on all DFs. An individual was reassigned to the other morph if its multivariate phenotype was closer to the centroid of that morph than to the centroid of its original morph.
Because of gaps in the data set, the density and length of glandular hairs on calyces and bracts were analysed separately with analysis of variance (ANOVA), based on the same model as was used in the MANOVA, with all factors treated as fixed.
Categorical data (plant colour, position of glandular hairs on bracts) were subjected to ordinal linear regression (position of glandular hairs on bracts), based on the same model as described above and a statistical testing procedure based on log-likelihood ratios. Analysis of plant colour was only performed for the common morph, since almost all fen plants fell into the same colour category (see Results).

Leaf collection and DNA extraction
Samples for molecular analyses were collected from the same fen morph populations on Gotland as those sampled for the common-garden study (SO1-SO3), whereas the common morph was represented by four populations on Gotland (A8, G2, G10, SA1), all from the same range as the fen populations, two populations on the Swedish mainland (A4, A7) and three populations on Saaremaa (EA1-EA3). Microsatellite data from seven previously studied populations of R. osiliensis (Rosi2-10) obtained from Talve et al. (2014a, b), were also included for reference purposes, giving a total of 19 sampled populations and molecular data from 509 individual plants (Fig. 1, Table 1). Leaves for the present study were collected from 40 plants per site, growing > 2 m apart. Leaves were dried and stored in plastic bags with silica gel until DNA extraction. Total genomic DNA was extracted using a modified CTAB-protocol (Doyle and Doyle 1990).

Molecular methods
Samples were prepared for microsatellite typing at 11 loci described by Ducarme et al. (2008), Houston and Wolff (2009)  PCR amplifications were carried out in a PTC-100 thermal cycler using the following program: initial denaturation at 94 °C for 5 min followed by 30 cycles of: 94 °C for 30 s, the annealing temperature (T a , Table 2) for 45 s and 72 °C for 45 s, followed by 9 cycles of: 94 °C for 30 s, 53 °C for 45 s and 72 °C for 45 s, and a final extension of 72 °C for 10 min (Jonstrup et al. 2018). PCR products were sent to Uppsala Genome Centre, where an internal size standard (Genescan, Liz500) was added and fragments were separated using an ABI3730XL DNA Analyser. Fragments were sized using the software GENIOUS 8.0 (Kearse et al. 2012) and quality-checked by eye. To enable exact comparison of our microsatellite data with those previously published for R. osiliensis (Talve et al. 2014a, b), eight of the R. osiliensis samples already analysed in Talve et al. were rerun together with the samples from Gotland. As it was still difficult to match the allele sizes at three loci (Ra58, Ra77 and Ra84), these loci were excluded from analyses including both Swedish and Estonian material.

Genetic structure and differentiation
To test whether the fen morph and the common morph on Gotland constitute distinct units and whether a genetic structure exists within the fen morph, we analysed the microsatellite data with the software STRU CTU RE 2.3.4, which identifies clusters that minimize Hardy-Weinberg and linkage disequilibria based on multilocus genotypes (Pritchard et al. 2000). A second analysis was performed to assess the relationship between the Gotland fen morph, the common morph of R. angustifolius from both Sweden and Estonia, and R. osiliensis from Saaremaa. We used settings according to recommendations given in the STRU CTU RE manual (Pritchard et al. 2000). Simulations were run for 1-10 clusters (K = 1 to K = 10) in the first analysis and for K = 1 to K = 24 in the second analysis. For each value of K, the analysis was repeated 20 times with a burn-in period of 100 000 iterations followed by additional 100 000 iterations to obtain parameter estimates. We ran the analyses with prior information on population identity and applied a model that assumes admixture between populations and correlated allele frequencies (Falush et al. 2003). The optimal K value, i.e. the number of clusters that best fit the data, for each analysis was determined by the Evanno method (Evanno et al. 2005), implemented in the web platform STRU CTU RE HARVESTER 0.6.94 (Earl Dent and von Holdt 2012). Separate runs for each K were aligned using CLUMPP 1.1.2 (Jakobsson and Rosenberg 2007), and the results were visualized with DISTRUCT 1.1 (Rosenberg 2004).
Patterns of differentiation between and within morphs were also evaluated using analysis of molecular variance (AMOVA), as implemented in ARLEQUIN 3.5.2.2 (Excoffier and Lischer 2010). We grouped populations by morph (the two different morphs of R. angustifolius with or without R. osiliensis as the third group) and ran two sets of analyses, based on pairwise F ST and R ST -values as measures of differentiation, respectively. F ST , also called the fixation index, is based on the infinite alleles model (IAM; Kimura and Crow 1964) and treats all distinct alleles as equally different, while R ST assumes a stepwise mutation model (SMM;Slatkin 1995;Rousset 1996) by taking into account the length of each microsatellite allele. F ST -based analyses suffer less from sampling variance than analyses based on R ST (Balloux and Goudet 2002). We visualized the genetic structure by means of principal coordinates analyses (PCOs), performed in the program NTSYS-pc 2.20 (Rohlf 1994) using pairwise F ST -values obtained from ARLEQUIN. To test for geographical structure at the regional scale, we performed a Mantel permutation test on geographical distances and F ST -values for pairs of populations on Gotland, using ARLEQUIN with data obtained in SPAGeDi 1.5 (Hardy and Vekemans 2002).

Genetic diversity within populations
The genetic diversity of individual populations was estimated as the number of alleles (A), the number of effective alleles (A e ), the number of private alleles (A p ), the observed and expected heterozygosity (H o , H e ) and the coefficient of inbreeding (F IS ), all calculated using GenAlEx 6.502 (Peakall and Smouse 2012). Since sample sizes differed between populations (Table 1), we also estimated allelic richness adjusted for sample size (A adj ) by means of rarefaction, as implemented in SPAGeDi 1.5 (Hardy and Vekemans 2002). Input files for analyses in STRU CTU RE, ARLEQUIN and SPAGeDi were created in CREATE (Coombs et al. 2008).

Phenotypic variation
Flowering time and overall plant morphology differed significantly between the fen morph and the common morph, and between the three populations representing each morph in the common-garden experiment (morph: F 4,414 = 81.52, P < 0.001; population: F 16,1265 = 39.36, P < 0.001; MANOVA). Judging from the DF plot, individuals of the fen morph generally had larger values along the first axis of separation (DF1, Fig. 2a, b) owing to the presence of more nodes, more branches and a later flowering start (Table 4). Environmental factors exerted significant effects on the mean phenotype (host: F 8,828 = 2.79, P < 0.01; block: F 12,1096 = 3.50, P < 0.001); however, there was no clearly visible pattern of response across the three host species along the two major axes of separation (Fig. 2b). Our analysis revealed a weak, though statistically significant, populationby-host interaction (F 32,1528 = 1.63, P < 0.05) but no morphby-host interaction (F 8,828 = 0.89, P > 0.05).
In the DFA of data pooled over blocks, host environments and replicate populations, 41% of the variation in flowering time and plant morphology could be attributed to differences between morphs (1-Wilks lambda = 0.41, P < 0.001). 81% of the individuals were correctly assigned to their original morph when they were re-classified according to the DFs (data not shown).
None of the individuals of either the fen morph or the common morph had normal or glandular hairs on the main surface of the calyx or corolla, but a thin zone of glandules and non-glandular hairs was always present along the margin of the calyx (A. Jonstrup, personal observation). Most of the glandules were sessile, but a few stalked glandular hairs were often present. Plants of the fen morph and the common morph differed in both the density and length of glandular hairs on the calyx (P < 0.001-0.01, ANOVA; Table 5) with fen morph plants having a denser cover (mean ± SD 10.8 ± 8.7 vs. 5.7 ± 5.1 for the common morph) and longer hairs (3.6 ± 3.1 vs. 3.0 ± 1.2 for the common morph); note, however, that there was considerable overlap between sd's in both cases. Populations within morphs differed significantly (P < 0.001) and both characters were significantly influenced by host species (P < 0.05 in both cases, ANOVA; Table 5). Effect sizes for block and interaction terms were too weak to be declared significant (P > 0.05, ANOVA; Table 5).
Plants of the fen morph had glandular hairs over a larger proportion of the bracts (with 32% having more than half of the bract covered by hairs) than plants of the common morph (no plant with the entire bract covered by hairs; P < 0.001, ordinal logistic regression), whereas the effect of population within morph was nonsignificant (P > 0.05, ordinal logistic regression). Growth on different host species influenced the area covered by glandular hairs (P < 0.01, ordinal logistic regression), with individuals grown on Trifolium having the highest scores (data not shown). The length of the glandular hairs on the bracts was not significantly affected by any of the factors considered (P > 0.05, ANOVA; Table 5).
The fen morph was distinctively uniform in colouration, with 99% of the plants being red and 1% being green. Plants of the common morph were more variable (75% green, 25% red), with no consistent differences between blocks, populations and host environments (P > 0.10, ordinal logistic regression).

Molecular-genetic differentiation
We performed separate analyses in computer program STRU CTU RE (Pritchard et al. 2000) to explore the genetic structuring of populations from Gotland only, and to explore the genetic structuring of the combined material from Gotland, mainland Sweden and Estonia. The resulting patterns are illustrated graphically in Figs. 3 and 4, respectively, and we have chosen to present the results given by K = 2 to K = 4, since in each case these are the most relevant partitions relating to the hypothesis formulated in the introduction.
In the STRU CTU RE-analysis restricted to the fen morph and the common morph on Gotland (Fig. 3), the optimal partitioning was obtained at K = 4. Individuals from fen morph populations split into two distinct clusters, one comprising individuals from two southern populations (SO2, SO3) and another comprising individuals from a northern population (SO1), together with those from a nearby roadside population (SA1). Individuals from two southern common morph populations (A8, G2) formed a third cluster, and individuals from a northern common morph population (G10) were joined in a fourth cluster. For K = 3, individuals from fen morph populations remained separated into two clusters. Population SO1 was still grouped with the nearby common morph population SA1, while individuals from the remaining common morph populations formed their own cluster (Fig. 3). For K = 2, individuals from the northern fen morph population (SO1) continued to form a distinct cluster together with individuals from the SA1 population, while individuals from the northern common morph population (G10) constituted another distinct cluster. The remaining fen morph and common morph plants were about equally likely to belong to either of the two clusters (Fig. 3).
In the STRU CTU RE-analysis based on the whole material from Sweden and Estonia (Fig. 4), the optimal partitioning was obtained at K = 2. The two clusters corresponded to individuals from Sweden and Estonia, with exception of individuals from the Estonian population EA3, which clustered together with the Swedish individuals. Swedish and Estonian individuals of the common morph had mixed probability of belonging to any of the two clusters, while individuals of the Gotland fen morph and R. osiliensis were assigned to different regional clusters. For K = 3, the Swedish individuals separated into two clusters, with individuals from the two southern fen morph populations (SO2, SO3) forming one cluster and the remaining ones, including the northern fen morph population SO1, forming an admixture of the two Swedish clusters. Most of Estonian individuals were kept as a single cluster, whereas those of R. angustifolius, especially population EA3, continued to show some similarity to those from Sweden.
The PCO on F ST -values for population pairs on Gotland (Fig. 5) confirmed the pattern found in the STRU CTU REanalysis. Populations were spread out in the PCO plot (Fig. 5), more or less forming three clusters, one with the two southern fen morph populations (SO2, SO3), one with the northern fen morph population (SO1) and the nearby roadside population (SA1), and one cluster with the remaining common morph populations (A8, G2, G10). No significant correlation was found between geographical and genetic distance between populations on Gotland (r = 0.23, P > 0.05, Mantel test).
Similarly, PCO on data from all Swedish and Estonian populations showed the same large-scale geographical pattern (Fig. 6) as that seen in the STRU CTU RE-analysis. Swedish and Estonian populations tended to fall into separate clusters (except for EA3), the most distinct being the Estonian R. osiliensis populations. Gotlandic fen morph populations continued to form two separate clusters, in this case with the northern cluster including not only populations SO1 and SA1, but also the Estonian R. angustifolius population EA3.
In the AMOVA based on fen morph and common morph populations from Gotland, a large portion of the diversity was due to differentiation between populations within morphs (F ST : 37.7%; R ST : 44.7%; Table 6), while only a minor portion was due to differentiation between morphs (F ST : 7.2%; R ST : 0.5%; Table 6). In the AMOVA based on all Swedish and Estonian populations, the between-morph component of differentiation was higher (F ST : 18.2%; R ST : 13.3%; Table 6), but much of the variation was still explained by differentiation between populations within morphs (F ST : 24.7%; R ST : 29.1%; Table 6). Fig. 3 The probability of each individual of Rhinanthus angustifolius from Gotland being member of a particular cluster identified by STRU CTU RE, based on 11 microsatellite loci. Individuals are arranged according to their population and region. Colours denote estimated membership fraction. K = 4 gave the best fit, but results are also shown for K = 2-3

Genetic diversity within populations
Analysis of diversity at all 11 microsatellite loci studied in the Gotland material showed fen morph populations (SO1-3) to have lower values than common morph populations (Table 7). In particular, measures of numbers of alleles per locus, richness and heterozygosity were usually the lowest for either of the two southern fen morph populations (SO2, SO3) and the highest for one of the southern common morph populations (G2, A8). Judging from a closer examination of allele frequencies (Supplementary data), the fen morph populations contained different subsets of alleles found in the common morph, especially at locus Rhos8 (fixed in SO1) and at loci Ra73, Ra75, Ra87 and Rhos8 (fixed in both SO2 and SO3).
Analyses based on the eight loci whose alleles could be matched between Swedish and Estonian populations confirmed the trend seen in the analysis of all 11 loci-relatively low diversity in Gotlandic fen populations-but also a tendency towards reduced genetic diversity in the Estonian R. osiliensis populations (Rosi2-10, Table 8).
The number of private alleles (A p ) was low for all populations, regardless of whether the analysis was based on all 11 loci (A p = 0-0.7, Table 7) or the restricted set of eight 'matching' loci (A p = 0-0.4, Table 8). Coefficients of inbreeding (F IS ) were always positive, with a trend towards higher values for common morph populations of R. angustifolius (Tables 7, 8).

Discussion
Although adaptation ultimately occurs at the level of individual populations, it is often useful to identify more broadly defined ecotypes, i.e. groups of populations connected to Fig. 4 The probability of each individual of Rhinanthus from Sweden and Estonia being member of a particular cluster identified by STRU CTU RE, based on a subset of 8 microsatellite loci. Individuals are arranged according to their population and region. Colours denote estimated membership fraction. K = 2 gave the best fit, but results are also shown for K = 3-4 particular habitat types (Turesson 1922), especially if the aim is to conserve the genetic diversity and long-term adaptive potential of species at the intraspecific level (Jonsell 2004;Lundqvist et al. 2007). Decisions on taxonomic status and conservation value are, however, only meaningful if ecotypes correspond to biologically meaningful entities at both the genetic and phenotypic level (Lowe et al. 2004;Mace 2004;Stuessy et al. 2014).
In this study, we performed a common-garden experiment and obtained microsatellite data to evaluate the distinctness and taxonomic status of populations of R. angustifolius found in calcareous spring fens on the Baltic island of Gotland. We compared them with other conspecific populations (collectively called the common morph), and also with R. osiliensis, an endemic taxon known from calcareous spring fens on the island of Saaremaa. Our experimental results are consistent with the presence of a phenotypically distinct fen ecotype (sensu Turesson 1922), characterized by red colouration as well as a later flowering start, a larger node and branch number, and a denser cover of glandular hairs on bracts and calyces than normally found in R. angustifolius. The Gotland fen ecotype is genetically different from R. osiliensis but appears to be closely related to conspecific populations of the common morph from Gotland. Based upon these and other observations, we propose a varietal status for the Gotlandic fen ecotype and give some recommendations for the conservation of this taxon.

Phenotypic differentiation
In accordance with previous observations (Lindell 2006), we found significant differences between the two morphs in flowering time and general plant morphology, similar to what has been observed in other rhinanthoid taxa (e.g. ter Borg 1972;Karlsson 1982;Zopfi 1993aZopfi , b, 1995. In particular, fen morph plants started flowering later and had a larger node number than individuals of the common morph, a direct consequence of the close developmental link between flowering start and node number in annual rhinanthoids (Wesselingh 2016). After pooling of data across environments, 41% of the total variation in phenology and gross morphology could be attributed to differences between the fen morph and the common morph; further, a majority of the plants (81%) were correctly re-classified to their original morph in the multivariate analysis. Thus, there was relatively little overlap between morphs across the host environments used in the common-garden experiment.
The percentage of correctly classified individuals was higher than the 61-73% correctly classified individuals in similar comparisons of two putative seasonal ecotypes within the common morph (ssp. angustifolius vs. ssp.    environments, two factors that may have resulted in a larger phenotypic overlap between populations in their experiment. In particular, in the present study, we did not use the stronggrowing legume Trifolium pratense as host species, a species that exerted strong and highly specific phenotypic effects on the plants used in the previous study (Jonstrup et al. 2016).
Consistent with work on other Rhinanthus species (de Soó 1929;de Soó and Webb 1972), we also identified differences in characters unrelated to flowering time and node number. First, almost all fen plants were distinctly red because of anthocyanin pigmentation, while a majority of the common morph plants had normal green colouration, irrespective of population and environment. Second, plants of the fen morph had a denser and more extensive cover of glandular hairs on the calyx and the bracts subtending the flowers, and the glandular hairs on the calyx margin were longer. There was, however, considerable overlap between morphs and the glandular hairs on the fen plants were never as long and dense as observed in R. osiliensis (A. Jonstrup, personal observation).
Calcareous spring fens on Gotland are characterized by cold upwelling water all year around and by a soil poor in available nutrients (especially phosphorous) due to high calcium content (Sjörs 1956;Pettersson 1958;Marklund 1982), two factors that could enhance the selective advantage of features characterizing the fen morph phenotype. For example, the high content of anthocyanin may act as a defence against oxidative stress caused by mineral imbalance (Landi et al. 2015). A similar hypothesis applies to the relatively dense cover of glandular hairs observed in the Gotland fen morph and, to an even greater extent, the densely pubescent R. osiliensis on Saaremaa. Glandular hairs of Rhinanthus have recently been shown to function as water excreting hydathode trichomes, which help to create a low water potential and thus to facilitate the transport of water and nutrients from the roots of the host plant (Světlíková et al. 2015). However, as the observed differences only concern structures directly associated with the flowers (calyces, bracts), more detailed studies are required to determine the function and adaptive significance of glandular hairs in the Gotlandic fen morph.

Differentiation at molecular marker loci
In contrast to the phenotypic differences found between the two morphs of R. angustifolius, there was no corresponding differentiation at the microsatellite loci investigated: the fen morph populations on Gotland were more strongly differentiated from one another than from geographically adjacent populations of the common morph, a pattern seen in both the PCO and STRU CTU RE analyses. Thus, only a modest part of the microsatellite variation (< 20%) could be attributed to differences between morphs, irrespective of whether the AMOVA was based on all 11 microsatellite loci studied in the Gotland material or the subset of eight loci studied in both Swedish and Estonian populations. These results are in agreement with other studies that found no, or only weak, molecular support for the delimitation of infraspecific taxa within Rhinanthus species (Houston and Wolff 2012;Pleines et al. 2013;Jonstrup et al. 2018).
Incongruence between phenotypic and molecular-genetic variation has several possible explanations. Houston and Wolff (2012) suggested that such incongruence found in R. minor could be a result of phenotypic plasticity. However, plasticity is not likely to be a major reason for the incongruence observed in the present study: environmental effects and their interactions with morph and population effects had minor influence on the mean phenotype in the commongarden experiment. Instead, it appears that the Gotlandic fen morph evolved in parallel from different populations of the common morph as an adaptive response to the same selection pressures, similar to what has been inferred for other plant species (Levin 2001) including rhinanthoids in other habitats (e.g. Euphrasia: Karlsson 1974Karlsson , 1976. Consistent with this scenario, the northern and southern fen populations contained different subsets of alleles found in the common morph, and more specifically, there was a close genetic similarity between the fen population SO1 and the common morph population SA1. Because we performed a common-garden experiment, we infer that the differentiation between the fen morph and the common morph of Rhinanthus has a heritable basis, and that they accordingly constitute true ecotypes. However, our data provide no evidence as to whether the observed difference has a pure genetic background or whether epigenetic modifications of the genome also contribute. Epigenetic modifications have been shown to influence phenotypic traits in many plant systems, and these modifications can be stable over multiple generations (Pikaard and Scheid 2014). For instance, epigenetic factors seem to explain much of the phenotypic variation expressed within genetically homogeneous accessions of the inbreeding annual Arabidopsis thaliana (Zhang et al. 2018). Epigenetic modifications can also contribute to morphological differentiation between geographical populations and, in some cases, taxonomically recognized entities (e.g. Paun et al. 2010). Nevertheless, studies linking ecotypes in the original sense of Turesson (1922) to epigenetic mechanisms seem to be rare.
The genetic divergence between populations was generally large (mean F ST ≈ 0.3; Online Resources 1-2), as commonly observed for annual plants with gravity-dispersed seeds (mean G ST = 0.38; Hamrick and Godt 1996), indicating limited gene flow between populations. However, a history of gene flow is still possible between populations in close proximity, such as the two southern fen populations on Gotland (SO2, SO3) and the aforementioned SO1 and SA1 populations. Regarding SO2 and SO3, which shared the same subset of alleles found in the common morph, we cannot exclude a common origin, perhaps from a more ancient lineage of R. angustifolius, as the cause of this and other similarities. Furthermore, we note that common morph populations from both Sweden and Estonia were separated by relatively short genetic distances. Most of these populations represent managed grassland sites (pastures, road verges and hay fields) and they may accordingly be connected through seed dispersal by machinery, animals or hay transport, both within and between countries (ter Borg 1972). The most striking case concerns the Saaremaa population EA3, which clustered with the Swedish populations (especially those from Gotland) rather than with the other Estonian populations (Fig. 6). This connection could reflect fairly recent transport by humans, for example in hay with cattle, as previously invoked to explain the genetic proximity between Icelandic and Danish populations of R. minor (Vrancken et al. 2012).

Genetic diversity and conservation biology
Calcareous spring fens on Gotland are unique environments with a characteristic flora that includes several rare and redlisted species, including cold-adapted taxa with their only Swedish occurrences in this habitat, e.g. Tofieldia calyculata and the local endemic Euphrasia salisburgensis var. schoenicola (Sjörs 1956;Pettersson 1958). Many fens on Gotland have been lost or have become increasingly fragmented as results of drainage and exploitation during the nineteenth and twentieth centuries (Marklund 1982;Martinsson 1997). General legislations protect wetlands in southern Sweden from further exploitation, but increasing nutrient supply due to atmospheric nitrogen deposition and overgrowth still pose a threat to spring fen vegetation (Martinsson 1997). These changes call for further actions to prevent habitat degradation and raise the question of how to best preserve the genetic diversity and evolutionary potential of taxa in this habitat (Moritz 1994;Frankham 2010).
The present study suggests that the Gotland fen morph of R. angustifolius constitutes a distinct spring fen ecotype that evolved locally on Gotland and that the fen morph plants fall into two genetic clusters, one containing two southern populations and another containing a northern population (together with a nearby common morph population). Thus, to conserve the genetic diversity of R. angustifolius in general, and the fen ecotype in particular, it seems reasonable to treat the southern and northern fen populations as separate management units, each with a high conservation value (Moritz 1994;Frankham 2010). In this regard, we also call for further studies of the SA1 population to assess the taxonomic status and conservation value of this population.
Gotlandic fen populations were not particularly inbred but had lower genetic diversity and fewer private alleles than other populations, especially those of the common morph. Low genetic diversity could reduce fitness and the longterm evolutionary potential of populations by increasing the expression of recessive deleterious alleles and by decreasing the heritable variation available to selection (Young et al. 1996;Storfer 1996). Microsatellite markers are, however, assumed to be selectively neutral with little or no phenotypic effect, making it difficult to draw firm conclusions about population viability (Frankham 2010). Further studies of the fen ecotype are needed to determine whether the low genetic diversity is associated with any reduction in population viability (Leimu et al. 2006) and whether present-day populations of the fen morph have been purged of alleles with negative effects on survival and fecundity (Crnokrak and Barrett 2002).

Taxonomic treatment
Although fen morph populations on Gotland are provisionally treated as R. osiliensis by Swedish conservation authorities (ArtDatabanken 2015), there is a clear separation between fen populations on Gotland and Saaremaa (Talve et al. 2014a, b). Our current results, based on microsatellite data from both Swedish and Estonian populations, support the hypothesis of Talve et al. (2014a, b) that fenadapted populations on Gotland and Saaremaa have separate origins-extending the parallel-adaptation scenario proposed above for individual populations of the fen morph on Gotland. These considerations, together with the patterns observed at the phenotypic level, suggest that the Gotlandic fen morph warrants recognition as a separate taxon. Since the fen morph contains a subset of the variation found within the common morph on Gotland, not only at putatively neutral marker loci but also in phenotypic characters such as plant colour, there should be great potential for rapid and recurrent shifts towards the fen morph, as a result of selection operating on pre-existing genetic variation within ancestral populations (Barrett and Schluter 2008;Pleines et al. 2013;Jonstrup et al. 2016). Thus, it seems reasonable to distinguish the Gotlandic fen morph as an infraspecific taxon within R. angustifolius rather than as a separate species. On the basis of guidelines given in the ongoing Flora Nordica project Jonsell (2004), varietal status is deemed most appropriate on the ground that it originated independently on at least two occasions on Gotland and that the different populations are most closely related to geographically adjacent populations of the common grassland morph rather than to each other (Hamilton and Reichard 1992;Jonsell 2004). Other floras often restrict the usage of infraspecific categories to the level of subspecies (e.g. de Soó and Webb