Plant intraspecific functional trait variation is related to within‐habitat heterogeneity and genetic diversity in Trifolium montanum L.

Abstract Intraspecific trait variation (ITV), based on available genetic diversity, is one of the major means plant populations can respond to environmental variability. The study of functional trait variation and diversity has become popular in ecological research, for example, as a proxy for plant performance influencing fitness. Up to now, it is unclear which aspects of intraspecific functional trait variation (iFDCV) can be attributed to the environment or genetics under natural conditions. Here, we examined 260 individuals from 13 locations of the rare (semi‐)dry calcareous grassland species Trifolium montanum L. in terms of iFDCV, within‐habitat heterogeneity, and genetic diversity. The iFDCV was assessed by measuring functional traits (releasing height, biomass, leaf area, specific leaf area, leaf dry matter content, Fv/Fm, performance index, stomatal pore surface, and stomatal pore area index). Abiotic within‐habitat heterogeneity was derived from altitude, slope exposure, slope, leaf area index, soil depth, and further soil factors. Based on microsatellites, we calculated expected heterozygosity (He) because it best‐explained, among other indices, iFDCV. We performed multiple linear regression models quantifying relationships among iFDCV, abiotic within‐habitat heterogeneity and genetic diversity, and also between separate functional traits and abiotic within‐habitat heterogeneity or genetic diversity. We found that abiotic within‐habitat heterogeneity influenced iFDCV twice as strong compared to genetic diversity. Both aspects together explained 77% of variation in iFDCV (Radj2 = .77, F 2, 10 = 21.66, p < .001). The majority of functional traits (releasing height, biomass, specific leaf area, leaf dry matter content, Fv/Fm, and performance index) were related to abiotic habitat conditions indicating responses to environmental heterogeneity. In contrast, only morphology‐related functional traits (releasing height, biomass, and leaf area) were related to genetics. Our results suggest that both within‐habitat heterogeneity and genetic diversity affect iFDCV and are thus crucial to consider when aiming to understand or predict changes of plant species performance under changing environmental conditions.

Intraspecific trait variation (ITV) depends on the available phenotypic trait plasticity of individuals within a population.
Phenotypic plasticity, that is, the phenotypic variation expressed by a single genotype under different environmental conditions (Hufford & Mazer, 2003;Nicotra et al., 2010;Sultan, 2000), might be one of the most important mechanisms for plants in reacting to environmental changes (e.g., land use, climate change; Agrawal, 2001;Arnold, Kruuk, & Nicotra, 2019;Gratani, 2014;Via et al., 1995). In general, environment and genetics can generate ITV (de Bello et al., 2011;Violle et al., 2012). The complex relationships between population-based ITV of functional traits and environmental heterogeneity of the habitats where traits of populations have been investigated, however, have not yet received much attention.
Considerable functional differences may provide improved resource partitioning due to differences in niche exploitation and/or a more flexible response to environmental changes (see Bucher et al., 2016;MacArthur & Wilson, 1967;Schweiger et al., 2018;Simpson, 1949;Violle et al., 2012). Within populations, increased ITV based on within-habitat heterogeneity should be able to increase adaptability with positive consequences for growth, reproduction, and survival. Environmental heterogeneity within habitats may thus lead to an increased number of different functional phenotypes and thus enhances ITV. Moreover, habitat heterogeneity is expected to influence the genotypic range of variation within a habitat: Variable environments can exert different selective pressures generating genetic heterogeneity (Gratani, 2014;Linhardt & Grant, 1996;Sakaguchi et al., 2019). Nevertheless, within a habitat of a population, environmental differences are usually lower and gene flow more frequently (for example due to missing geographical barriers) than between habitats across larger scales. Within-habitat heterogeneity might also enhance the occurrence of different genotypes due to different resource exploitation possibilities, increasing genetic variation (see, e.g., Agashe & Bolnick, 2010;Reusch, Ehlers, Hämmerli, & Worm, 2005). Therefore, within-habitat heterogeneity affects ITV directly and genetic diversity indirectly. However, ITV, within-habitat heterogeneity, and genetic diversity may interact in complex ways under natural environmental conditions. For example, plasticity of traits (generating ITV) is able to influence the selective effect of within-habitat heterogeneity on genetic diversity, whereas connectivity and dispersal (gene flow) among habitats can affect selection on genetic diversity (Ghalambor, McKay, Carroll, & Reznick, 2007;Linhardt & Grant, 1996;Reisch & Schmid, 2019;Vellend & Geber, 2005).
In this study, we aim to investigate the relative effects of abiotic within-habitat heterogeneity and genetic diversity on intraspecific trait variation (ITV). Investigations are based on 260 individuals from 13 Central European populations of the rare (semi-)dry calcareous grassland species Trifolium montanum L. (mountain clover; Figure 1). We addressed the following questions: Is ITV related to abiotic within-habitat heterogeneity and/or genetic diversity? If so, to what extent is ITV explained by either aspect? Which functional traits are related to abiotic within-habitat heterogeneity and/or genetic diversity?  (Jäger, 2011), but they also occur along with shrub and forest margins. In Central Europe, the species is quite rare because of degradation and fragmentation of (semi-)dry grasslands (Garve, 2004;Schleuning & Matthies, 2008;Schleuning et al., 2009). Trifolium montanum is diploid with 2n = 16 (Rice et al., 2015).

| Study locations and sampling
We focused on 13 locations in total, from which 12 are situated in Germany and, to cover a larger range of environmental conditions, one in Austria (Table 1). Trifolium montanum populations at the 13 different locations are independent of each other: The average distance of mating events in this species is quite low (10 m), and thus, the majority of mating events occur on small distances (pollen dispersal up to a distance of 324 m possible; Matter, Kettle, Ghazoul, Hahn, & Pluess, 2013). This indicates a reduced potential for pollen-mediated long-distance dispersal. For two Trifolium species, the estimated long-distance dispersal via seeds is also estimated to be only six to 10 m (Vittoz & Engler, 2007). We estimated an average distance among study locations of ca. 133 km (80 km without the distant location KW) and standard deviation of 124 km (41 km without the distant location KW, see Table S3 for details and the applied "geosphere" r package vers. 1.5-5 (Hijmans, 2016) using the Vincenty ellipsoid method). Distances are thus too large for direct pollen exchange or seed dispersal. There is only the possibility for direct gene flow between locations Bo and Ha, and Ba and St. However, locations are separated by large agrarian areas (particularly Ha), forests, and roads, minimizing the probability of pollen exchange or seed-mediated long-distance dispersal and thus the potential of gene flow. Trifolium montanum is consumed by grazing animals, but Bo was grazed by sheep and goats, St by cattle, and the management of Ha is unknown, suggesting rather very local grazer movements instead of habitat connection via transhumance.
From each location, we collected 20 individuals of a T. montanum population, totaling 260 individuals. We attempted to distribute individual sampling points equally within a habitat. Fieldwork was done in July 2015 (Austria) and from May to June 2016 (Germany).

| Functional traits-measurements and ecological meaning
As traits change with the season Römermann, Bucher, Hahn, & Bernhardt-Römermann, 2016), we started sampling in lowlands and finished sampling in higher altitudes (compare also Tautenhahn, Grün-Wenzel, Jung, Higgins, & Römermann, 2019). To ensure the comparability of functional traits among populations, we only sampled flowering and early fruiting individuals to control for phenology . All functional traits were measured on 20 different individuals per population, and all leaf functional traits were measured on two leaves per individual.
In the field, we measured F v /F m and PI on absorption basis after 30 min by high intensity focused LED (3,500 µmol/m 2 * s −1 intensity and wavelength peak 627 nm) with Pocket PEA, preparing leaves according to the manufacturer's instructions (Hansatech Instruments Ltd., King´s Lynn, England;Strasser, Srivastava, & Tsimilli-Michael, 2000;Strasser, Tsimilli-Michael, & Srivastava, 2004). The ratio of variable fluorescence to maximal fluorescence (F v /F m ) is related to the efficiency of PS II electron transport and indicates abiotic and/or biotic stress due to photoinhibition (Butler & Kitajima, 1975 1976). The performance index (PI) represents the photosynthetic performance of a chlorophyll molecule, the vitality of the plant, and its ability to resist constraints from outside (Bucher, Bernhardt-Römermann, & Römermann, 2018;Clark, Landolt, Bucher, & Strasser, 2000;Strasser et al., 2000). We also determined releasing height (RH) as the shortest distance between the ground and the highest flower head [m] (Cornelissen et al., 2003), and cut total fresh aboveground biomass directly above taproot.
F I G U R E 2 (a) Distribution range of Trifolium montanum in Europe (light gray) according to Meusel and Jäger (1998). The black square indicates the sampling area in Central Germany. The black dot represents the sampling location in Austria ("KW"). (b) Sampling scheme of the present study in Central Germany (see Table 1  and its water-saturated fresh mass [g]. Specific leaf area (SLA) tends to be positively correlated with potential relative growth rate, but negatively with leaf dry matter content (LDMC) that represents leaf longevity/robustness (Cornelissen et al., 2003;Pérez-Harguindeguy et al., 2013;Römermann et al., 2016).
We took stomata imprints using nail polish from two leaves per individual. Utilizing an Olympus CH40 microscope, we counted stomata density at 200× magnification and measured guard cell length and width at 400× magnification. We determined individual mean values of four measurements for stomata density and of eight measurements for stomata length and width measurements of the abaxial leaf surface. We also calculated stomatal pore surface (SPS) as guard cell length [µm] * guard cell width [µm * π * 4 -1 (Balasooriya et al., 2009) and stomatal pore area index (SPI) as the product of (guard cell length) 2 [mm 2 ] * stomatal density [1/mm 2 ] (Sack, Cowan, Jaikumar, & Holbrook, 2003). SPS characterizes stomata size, SPI indicates stomatal conductance, and both traits are known to change along abiotic environmental gradients (Bucher et al., , 2017Woodward, Lake, & Quick, 2002).

| Assessment of habitat characteristics
We characterized each location with a maximum of five environmental replicate measurements. The replicates (grid cell 2 m × 2 m) were equally distributed within the local range of each population (20 to 4,600 m 2 , unpubl. data). Therefore, the distance between replicates differed according to habitat sizes. For one location (Eh, see Tables   1 and 2), we reduced the number of replicates to four due to the limited habitat size (~20 m 2 ). In total, we analyzed n = 64 records.

| Population genetics and laboratory work
We isolated DNA of sampled individuals from approximately 20-25 mg dry leaf material using a modified CTAB protocol (Doyle TA B L E 1 Location, date of sampling, latitude (lat. (N), longitude (long. (E.), and mean coefficients (with standard errors in brackets) for variation of intraspecific functional trait variation (iFD CV ), abiotic within-habitat heterogeneity (HD), and mean genetic diversity (GD; H e ) of 13 Trifolium montanum populations TA B L E 2 Coefficients of variation (CV) of particular functional traits (n = 260 individuals), abiotic factors (n = 64 replicates), and population genetic indices (n = 255 individuals) based on nine microsatellite markers (Matter et al., 2012)

| Data analyses
All statistical analyses (except genetic diversity calculations) were performed with R vers. 3.6.0 (R Core Team, 2019). We calculated means for numerical variables, medians for the ordinal variable slope exposure, and coefficients of variation (CV) for diversity variables as the ratio of standard deviation to mean. To evaluate the multiple linear regression model, we exceptionally used the adjusted coefficient of determination (R 2 adj ) instead of the coefficient of variation (R 2 ). This has the benefit of avoiding model overfitting in R 2 calculation (see Crawley, 2007).

| Intraspecific functional trait variation (iFD CV )
To detect erroneous entries (errors in measurement) in functional traits, we excluded all records (seven trait measurements in total) from the dataset with a distance of >4 standard deviations from the mean of all individuals (compare Díaz et al., 2016;Kattge et al., 2011). We deleted F v /F m outliers and the respective PI values (identical source, Pocket PEA). We also checked collinearity among traits, that is, when two or more traits were highly correlated (r > ~│.7│; Dormann et al., 2013). We assessed correlations among functional traits with Spearman's rank coefficient (r SP , cor.test()) due to non-normal distribution of data.
ITV (or functional diversity, "FD") can be calculated in different manners. Functional trait dissimilarity among or within species is used to calculate FD via trait distance matrices and dendrograms (Petchey & Gaston, 2006;Tilman, 2001; see in-

| Within-habitat heterogeneity (HD)
We tested for correlations among environmental factors with Spearman's rank coefficient (r SP , cor.test()) due to non-normal distribution of data. The r package "corrgram" vers. 1.13 (Wright, 2018) was used to visualize correlations. As explained above for functional traits, we checked for collinearity (r > ~│.7│; Dormann et al., 2013) and excluded the factors CEC K , CEC Ca , pH KCl , C org , CaCO 3 , T a , P a (r SP > .7), and CEC Mg (r SP = .50), and one variable due to an almost complete absence of variation (CEC Na ; see Figure S2 for correlations). Afterward, we calculated abiotic within-habitat heterogeneity (HD) as location-wise mean CV of altitude (CV altitude ), slope exposure (CV slope exposure ), slope (CV slope ), leaf area index (CV LAI ), soil depth (CV soil depth ), potential soil cation-exchange capacity (CV CECpot ), pH (CV pH ), soil nitrogen content (CV N ), soil phosphor content (CV P ), and soil potassium content (CV K ; see Tables 1 and 2).

| Genetic diversity (GD)
We scored microsatellite fragments with an internal size standard. We proved the scoring procedure at least three times and removed ambiguous results. Analyses were conducted for all individuals characterized by at least four microsatellite loci resulting in a final sample size of n = 255 individuals (see Table 2). Mean loci coverage was 90% per individual, that is, in mean, 90% of loci were present in an individual. To ensure that sampled populations represent the same genetic line of the species, we calculated individual-and population-wise distance matrices based on Nei´s genetic distance (Nei, 1978) and conducted principal coordinate analyses (PCoAs) in GenAlEx vers. 6.503 (Peakall & Smouse, 2006. Moreover, we performed analyses in STRUCTURE vers. 2.3.4 (Pritchard, Stephens, & Donnelly, 2000) setting an admixture model (with correlated allele frequencies), burn-in to 5,000, MCMC to 50,000, and K to one to 13 (10 replicates per K). The optimal K was determined by STRUCTURE HARVESTER (Earl & vonHoldt, 2012) using the Evanno method. We merged the replicates of each optimal K (K = 2 and K = 9) with CLUMPP vers. Moreover, relationships among iFD CV , HD, and GD in KW fit those observed among Central German populations (see Table 1 and Population size of T. montanum was positively related to genetic diversity (Karbstein et al., unpubl. data), and varying genetic diversity is needed to examine its effect on iFD CV .

| Saturation of iFD CV , HD and GD
To examine whether iFD CV , HD, and GD (see Statistical modeling section) are saturated within populations/habitats, we used an r script (Karbstein, 2020) that randomly chose one to 20, one to five and one to 18 (18 individuals genotyped per population at a minimum) samples and calculated iFD CV , HD and GD, respectively with 100 iterations per step. Mean iFD CV , HD and GD values were plotted against sample size for each population/habitat. We used r functions implemented in the packages "dbplyr" vers.

| Statistical modeling
To assess whether iFD CV is related to within-habitat heterogeneity and genetic diversity (and to assess which specific genetic diversity index best explains iFD CV ), we employed a multiple linear regression model with iFD CV as the dependent variable, and HD and genetic  Table 2 for abbreviations. Significance levels: ***p < .001 and *p < .05 multiple linear regression model contained HD and expected heterozygosity (H e ) as explaining, independent variables (see Table S1).
H e is widely used as genetic diversity index, and it less depends on population history (e.g., bottlenecks) compared to the other indices (Freeland et al., 2011;Kalinowski, 2004;Szczecińska, Sramko, Wołosz, & Sawicki, 2016). Therefore, we used population-wise F I G U R E 4 Relationships between coefficient of variation of particular traits (CV traits ) and abiotic within-habitat heterogeneity (HD) in 13 Trifolium montanum populations (n = 260 individuals) of Central Europe. 95% confidence intervals are drawn for all (marginal) significant relationships. Dotted regression lines represent only marginal significant relationships (0.05 < p < .1). See Table S1 for detailed model statistics, and Table 2  To illustrate separate relationships between iFD CV and HD or GD (see Figure 3), we executed two linear regression models and plotted regression results. Additionally, we conducted a F I G U R E 5 Relationships between coefficient of variation of particular traits (CV traits ) and genetic diversity (GD, H e ) in 13 Trifolium montanum populations (n = 255 individuals) of Central Europe. 95%-confidence intervals are drawn for all (marginal) significant relationships. Dotted regression lines represent only marginal significant relationships (.05 < p < .1). See Table S1 for detailed model statistics, and Table 2 for abbreviations. Significance levels: *p < .05 and '=.1 > p> .05 linear regression model to assess whether HD is associated with GD. GD was handled as dependent and HD as the independent variable.
To assess which functional traits were related to HD and/or GD (see Table S1), we performed linear regression models with log-transformed trait variation as dependent variable and HD or GD as independent variable. We log-transformed the dependent model variable to achieve normality and/or linearity, and checked normality with the Shapiro-Wilcox test. Model assumptions were again visually examined as described above.
To test for spatial autocorrelation among populations/habitats (e.g., closer populations with more similar genetic diversity and more similar iFD CV ) within linear regression models, we calculated Moran's I (Moran, 1950) values using the R function correlog() function included in the r package "ncf" vers. To support the interpretation of results found between iFD CV and HD, we also examined correlations among particular functional traits and abiotic environmental factors chosen to calculate HD. We added a value of 1 to all CVs and logarithmized traits and environmental factor CVs (CVs have to be >1) to achieve normal distribution.
We used the rcorr() function within the r package "Hmisc" vers. 4.2-0 (Harrell, 2019) to calculate a correlation matrix based on Pearson's rank correlation coefficient. We carried out the corrplot() function implemented in the r package "corrplot" vers. 0.84 (Taiyun & Simko, 2017) to visualize the correlation matrix only considering p values below .1.

| Relationships of functional traits with environmental factors and genetic diversity
Most functional traits were significantly positively related to HD and some to GD (Figures 4 and 5; see Table S1 for detailed statistics). We observed a higher number of (marginally) significant relationships between functional traits (CV RH , CV AGB , CV SLA , CV LDMC , CV Fv/Fm, and CV PI ) and HD than between functional traits (CV RH , CV AGB and CV LA ) and GD. CV RH is strongest related to HD (R 2 = .73), followed by CV SLA (R 2 = .39), CV Fv/Fm (R 2 = .39), CV PI (R 2 = .39), CV AGB (R 2 = .37), and CV LDMC (R 2 = .28). With GD, CV LA (R 2 = .34), CV AGB (R 2 = .32), and CV RH (R 2 = .29) exhibited significant relationships of similar strength.
We found (marginally) significant positive correlations between functional traits and abiotic environmental factors chosen to calculate HD (Figure 6, see Table S2): Some trait CVs are positively correlated to CV slope exposure (CV RH , CV AGB , CV LA, and CV SPI ) and CV slope (CV RH , CV Fv/Fm and CV PI ). Many traits were positively correlated to CV of soil nutrients, that is, CV N (CV RH and CV AGB ), CV P (CV SLA and CV LDMC ), and CV K (CV Fv/Fm and CV PI ). No significant trait CV correlations were found to CV altitude , CV LAI , CV soil depth , CV CECpot, and CV pH.
Correlations between trait CVs and environmental CVs revealed widely positive coefficients (~71%).

| D ISCUSS I ON
Connecting intraspecific functional trait variation with genetic and environmental variation is an important ecological challenge.
This study showed that population-wise intraspecific functional trait variation (iFD CV ) of T. montanum can be attributed to a high extent (77%) to both abiotic within-habitat heterogeneity and population genetic diversity under natural environmental conditions ( Figure 3). Interestingly, within-habitat heterogeneity statistically affected iFD CV considerably stronger than genetic diversity (Figures 4 and 5).

| Relationships between intraspecific functional trait variation and within-habitat heterogeneity
Variation of morphology-related functional traits RH, AGB, and LA was mainly correlated to variation of habitat slope exposure and slope, whereas variation of (eco-)physiology-related traits (SLA, LDMC, F v /F m , PI, SPS, and PCI) was predominantly correlated to variation of slope characteristics and soil factors (Figure 6). Within a habitat, different slopes and slope exposures, influencing soil humidity, may have led to an increase or reduction in height, biomass, and leaf area of T. montanum. For example, drought stress is known to limit nutrient uptake and thus photosynthesis and plant growth (Farooq, Wahid, Fujita, & Basra, 2009;Jaleel et al., 2009;Pérez-Harguindeguy et al., 2013). Also, RH and AGB were positively associated with soil nitrogen content (N), and N-deficient soils in T. montanum habitats probably constrained height and biomass accumulation of individuals as well (see, e.g., Ågren, Wetterstedt, & Billberger, 2012). (Positive) associations between plant height or biomass to soil properties and particularly N are in line with literature (Cornelissen et al., 2003;Razaq, Zhang, & Shen, 2017;Reich & Hobbie, 2013). In contrast, variation of SLA and LDMC was correlated to soil phosphorous content (P). Phosphor regulates protein biosynthesis and development of new plant tissue (see Kerkhoff, Fagan, Elser, & Enquist, 2006), and P-deficient soils may thus impact relative growth rate and leaf robustness. SLA and LDMC respond also to soil properties (Cornelissen et al., 2003). P has been shown to form only rare and weak positive intraspecific relationships with SLA and particularly in herbaceous species when N is in abundant supply (wild rice; Dwyer, Hobbs, & Mayfield, 2014;Sims, Pastor, & Dewey, 2012). Moreover, soil potassium content (K), besides slope affecting soil humidity (see above), influenced performance and vitality (PI and F v /F m ). High soil potassium content was mainly found in CaCO 3 -rich soils, which T. montanum prefers (Jäger, 2011). Due to physiological constraints, low K and CaCO 3 conditions may decrease the performance/vitality of T. montanum individuals while high K and CaCO 3 may increase it. Studies already revealed that Intraspecific functional trait variation captured by iFD CV is probably the response of populations to increased abiotic and biotic environmental differences within their habitats (see also Ghalambor et al., 2007;Nicotra et al., 2010Nicotra et al., , 2015: as shown, iFD CV was positively correlated to within-habitat heterogeneity, suggesting that the more environmentally variable a habitat, the higher the intraspecific functional trait variation in T. montanum populations. Results are line with literature showing associations between functional trait values and diversity, and environmental conditions within and across species (e.g., Albert et al., 2010;Bernhardt-Römermann, Gray, et al., 2011;Bucher et al., 2016;Díaz & Cabido, 2001;Gratani, 2014;Karbstein et al., 2019;König et al., 2018;Violle et al., 2007). small-scale habitat differences). Thus, iFD CV likely affects population growth and reproduction with positive consequences for survival and fitness (see also Nock et al., 2016;Violle et al., 2007).

| Relationships between intraspecific functional trait variation and genetic diversity
Genetic diversity, in terms of microsatellite variation, was positively related to intraspecific functional trait variation (iFD CV ; Figures 3b and 5). This observation is also in line with literature indicating F I G U R E 6 Visualized correlation matrix based on Pearson correlation coefficients between variation of particular traits (CV trait ) and particular abiotic environmental factors (CV factor ) in 13 Trifolium montanum populations (n = 260 individuals) of Central Europe. We only illustrated (marginal) significant results (see Results). Width of an ellipse reflects the correlation coefficient, that is, the higher a correlation coefficient (in positive and negative direction), the narrower the ellipse. See Table 2 for abbreviations, and Table S2 for statistics. Significance levels: **=p < .01, *=p < .05 and '=0.1 > p > .05 positive relationships between trait variation/plant fitness (iFD CV , plant fitness, see above) and genetic diversity within species (e.g., Leimu et al., 2006;Waitt & Levin, 1998). Microsatellite markers are frequently applied to capture population genetic diversity (e.g., Matter et al., 2013;Matter et al., 2012;Prinz, Weising, & Hensen, 2009). They are widely distributed throughout genomes, and while regulatory functions in gene expression are known, these markers are presumed to predominantly occur in non-coding regions and thus to be under neutral selection (see also Ellegren, 2004;Li, Korol, Fahima, Beiles, & Nevo, 2002;Vieira, Santini, Diniz, & Munhoz, 2016). The weak relationship of iFD CV with genetic diversity may be explained by the neutrality of applied microsatellite markers in relation to selected functional traits. However, our intention was not to explain functional traits with particular microsatellites but to assess whether iFD CV and/ or particular trait variation coincide with genetic diversity.
It is likely that genetic diversity limits and influences the range of iFD CV . Mutation and recombination events create genetic variation and thus novel functional trait variation within a population. Natural selection probably acts on the genetic basis of trait variation, which in turn probably affects the range of iFD CV . After modeling, expected heterozygosity (H e ) best-explained iFD CV . H e is based on allelic structure, represents genotype and allele frequencies, and is less sensitive to population history (Freeland et al., 2011;Kalinowski, 2004;Szczecińska et al., 2016). Heterozygosity within individuals and populations probably influenced iFD CV because it enhances the reaction norm and adaptability and thus affects intraspecific trait variation (Boulding, 2008;Freeland et al., 2011;Reed & Frankham, 2003). Interestingly, only the variation of morphology-related traits was associated with genetic diversity in T. montanum (see also Waitt & Levin, 1998). In contrast, both variation of morphology-and (eco-)physiology-related traits was linked to within-habitat heterogeneity. (Eco-)physiology-related traits (gas exchange and photosynthesis) tend to have a higher heritability (trait variation due to genetic variation) compared to morphology-related traits (morphology and vegetative performance, Geber & Griffen, 2003), and should thus be more sensitive to genetic variation . However, an explanation might be that heterozygosity effects (see, e.g., Boulding, 2008;Freeland et al., 2011;Reed & Frankham, 2003) are stronger pronounced in morphology-related traits leading to similar trait variation and genetic variation based on microsatellites.
However, a positive relationship between genetic diversity and iFD CV may be strengthened by the self-incompatible nature of T. montanum (see, e.g., Leimu et al., 2006;Reed & Frankham, 2003;Schleuning et al., 2009). Observations between classical fitness parameters and genetic diversity of self-incompatible species are frequent and can be explained by pollinator limitation in fragmented populations with a low density of flowering individuals (Leimu et al., 2006;Schleuning et al., 2009). This process enhances the loss of genetic diversity in smaller populations, and it extends the range of genetic variation between small and big populations, probably also affecting the range of IFD CV .
Thus, relationships between genetic diversity and iFD CV , which influences plant fitness directly and indirectly (Nock et al., 2016;Violle et al., 2007), might be strengthened in T. montanum.

| Differentiated view on relationships among intraspecific functional trait variation, within-habitat heterogeneity, and genetic diversity
Relations between iFD CV , within-habitat heterogeneity and genetic diversity (Figure 7) are reported from several species and discussed F I G U R E 7 A conceptual model of relationships among intraspecific trait variation (functional diversity; iFD CV ), abiotic within-habitat heterogeneity (HD), and genetic diversity (GD) in T. montanum. HD influenced iFD CV twice as much as GD symbolized by different circle sizes and arrow strengths. Results and percents are extracted from the multiple linear regression model (R 2 adj = .77, F 2, 10 = 21.66, p < .001). HD can lead to (reversible) short-term responses, that is, to modification of functional trait expression (phenotypic modifications, variation). In contrast, GD is controlled by selection and is a prerequisite for adaptation through natural selection. iFD CV thus also depends on the available genetic variation within a population. Habitat heterogeneity and genetic diversity are not significantly related in this study (R = .10, F 1, 11 = 2.37, p = .15; dashed line) in literature cited above. Under natural environmental conditions, within-habitat heterogeneity and genetic diversity probably act on iFD CV in complex ways. For example, several factors may influence whether environmental variation within a habitat promotes the occurrence of different genotypes. Although a positive relationship between within-habitat heterogeneity and population genetic diversity was expected from literature (particularly for self-incompatible, outcrossing species, like T. montanum, with moderate to high population genetic diversity), we did not observe a significant effect.
Several reasons can explain this result. Large phenotypic plasticity of traits can enable individuals to inhabit different environmental niches potentially shielding them from natural selection (e.g., Ghalambor et al., 2007). Thus, there would have been no need to select for higher functional trait variation, weakening the relationship between within-habitat heterogeneity and population genetic diversity. Moreover, some close T. montanum populations are potentially connected (or were at least in the past) due to sheep and goat grazing (transhumance). Connectivity would have allowed for gene flow among them (Linhardt & Grant, 1996;Reisch & Schmid, 2019;Vellend & Geber, 2005) superseding (partially) local genotypes and altering population genetic diversity and adaptation to local environmental habitat conditions. However, some T. montanum populations in lowly/highly variable habitats are characterized by low/high genetic diversity (and low/high iFD CV ) indicating an association between environment and genetics (and trait variation) within habitats of a species (see also Gram & Sork, 2001;Huenneke, 1991;Linhardt & Grant, 1996). Varying resource exploitation of different genotypes may explain the observed pattern (see, e.g., Agashe & Bolnick, 2010;Reusch et al., 2005). Moreover, within-habitat heterogeneity could be underestimated in large T. montanum habitats (e.g., Er, If, Bo, St) characterized by large population sizes (Karbstein et al., in prep.), genetic diversity, and iFD CV . More environmental replicates would have potentially led to higher within-habitat heterogeneity estimates additionally strengthening positive relationships among iFD CV, within-habitat heterogeneity, and genetic diversity.
Saturation of diversity variables (iFD CV , HD, and GD) is an important feature. Unsaturated variables can bias relationships, and can lead to false results and conclusions. In T. montanum, at least five to 10 samples were sufficient to saturate iFD CV within populations (see also Bastias et al., 2017). Within-habitat heterogeneity was saturated by three to five samples per habitat, and genetic diversity reached the plateau at 10 to 15 samples per population, not biasing our regression models.
In addition, epigenetic processes, like DNA methylation or activation of transposable elements, in response to environmental variation can also alter phenotypic plasticity (Nicotra et al., 2015;Weinhold, 2006) and thus ITV potentially explaining a particular amount of unexplained variation in regression models. Our primary goal was not to separate the environmental component from the genetic one but to understand the relative importance of both environment and genetics on iFD CV , and particularly, how populations react under natural environmental conditions. To clearly separate the environment from the genetic impact on iFD CV , common garden experiments under controlled environmental conditions are necessary. Moreover, genomic data will provide more insights into genetic variation of populations, and in investigating relationships between functional trait variation and genetic variation. Connecting intraspecific functional trait, environmental, and genetic variation remains still an important ecological challenge.

| Impact on biodiversity research
Trait variation is probably of major importance to plants short-term adjustment on (rapid) environmental changes (see also Arnold et al., 2019;Gratani, 2014). Genetic diversity influences the range of trait plasticity and thus trait variation within a population, which can be advantageous for short-term responses (e.g., land use abandonment) by offering genetic variants that are fitter under novel environmental conditions. Moreover, in the long term (e.g., considering anthropogenic climate change), genetic diversity offers variation for natural selection to act and thus allows for adaptation to novel habitat conditions.
Environmental habitat aspects and population genetics should be considered in biodiversity research dealing with intraspecific functional trait variation at population, community, and ecosystem level. Consideration of these aspects can prevent bias and misinterpretation of trait variation analyses, for example, comparing trait variation between sites where trait differences cannot be attributed to environment or genetics. Habitat features are directly extractable from field measurements as shown in this study or potentially from databases with a high spatial resolution (e.g., WorldClim;Hijmans et al., 2005). If genetic features of a species cannot be investigated due to a lack of suitable markers, population size can be used as a cautious proxy of genetic diversity (see Leimu et al., 2006).
Our study demonstrates the potential of deriving intraspecific functional trait variation based on environmental and genetic aspects (or its proxies) and provides empirical evidence to encourage the incorporation of intraspecific functional trait variation into interspecific comparisons (see also Albert et al., 2010;de Bello et al., 2011;Violle et al., 2012). Directly measured species-specific intraspecific functional trait variation, but also values from databases provide the possibility for a better understanding of community and ecosystem responses to environmental changes and a more realistic estimation of ecosystem functioning.

ACK N OWLED G M ENTS
The Institute of Ecology and Evolution (Friedrich-Schiller-University, Jena, Germany) financially supported this research. We acknowledge the lower nature conservation authorities ("UNB") in Arnstadt, Bad Salzungen, Jena, Rudolstadt, and Sondershausen for permitting access to conservation areas, and the Thüringer Landesanstalt für Landwirtschaft und Ländlichen Raum (TLLLR) for preparing and conducting soil analyses. We also thank Janin Nauman, Annika Lambert,