Changes in global DNA methylation under climatic stress in two related grasses suggest a possible role of epigenetics in the ecological success of polyploids

Polyploidization drives the evolution of grasses and can result in epigenetic changes, which may have a role in the creation of new evolutionary lineages and ecological speciation. As such changes may be inherited, they can also influence adaptation to the environment. Populations from different regions and climates may also differ epigenetically; however, this phenomenon is poorly understood. The present study analyzes the effect of climatic stress on global DNA methylation based on a garden collection of two related mountain grasses (the narrow endemic diploid Festuca tatrae and the more widely distributed mixed-ploidy F. amethystina) with different geographic ranges and ecological niches. A lower level of DNA methylation was observed for F. tatrae, while a higher mean level was obtained for the diploid and tetraploid of F. amethystina; with the tetraploids having a higher level of global methylated DNA than the diploids. The weather conditions (especially insolation) measured 24 h prior to sampling appeared to have a closer relationship with global DNA methylation level than those observed seven days before sampling. Our findings suggest that the level of methylation during stress conditions (drought, high temperature and high insolation) may be significantly influenced by the ploidy level and bioclimatic provenance of specimens; however an important role may also be played by the intensity of stress conditions in a given year.

(1) What changes in global DNA methylation are observed in specimens of F. amethystina and F. tatrae kept in common garden conditions, grown under stress (drought, high temperature and high insolation) and under more favorable conditions, according to species and ploidy level? (2) Which climatic factors, soil humidity or phenotypic characteristics of specimens measured directly during the experiment are most closely related to DNA methylation level? (3) Which characteristics of the original provenance of the plants are most closely related to the DNA methylation levels in a common garden conditions?

Results
Global DNA methylation in stress and control conditions. A higher level of median DNA methylation was observed in all groups (species and cytotypes) during climatic stress (i.e. thermal, drought and insolation) (Fig. 1). Higher levels of DNA methylation were observed in 2019: the median values did not exceed 1.5% in 2018; however, they did not fall below 1.2% the next year ( Fig. 1). Of the two studied species, F. amethystina demonstrated higher mean global DNA methylation values in both 2018 and 2019, with the highest values being achieved by the tetraploid F. amethystina. In 2018, statistically significant differences were only observed between F. tatrae and tetraploids of F. amethystina (Fig. 1); in 2019, only slight differences were observed between cytotypes.
The variability of the epigenetic reaction of individual plants, i.e. the difference in the degree of DNA methylation between stress and control conditions, differed between species and cytotypes (Fig. 2). In both years, the greatest decrease in DNA methylation from stress to control conditions was noted for the tetraploid F. amethystina (Fig. 2). A smaller decrease was observed for the diploids of both species; however, no differences in median DNA methylation were demonstrated in diploid F. amethystina in 2018. Both diploid and tetraploid specimens of F. amethystina demonstrated visibly greater changes in DNA methylation between stress and control conditions in 2019 compared to 2018 (Fig. 2).
For all individuals, significant (p < 0.05) but moderate positive correlations in global DNA methylation levels were found between stress and control conditions within each measurement year (Pearson's coefficient: 0.34 for 2018 and 0.40 for 2019), i.e. individuals with a higher level of DNA methylation under stress conditions also maintained a higher level of methylation under control conditions. Consequently, the differences in methylation in a given year between stress and control were negatively correlated with the level of methylation during stress (− 0.85 for 2018 and − 0.46 for 2019): those specimens with higher methylation during stress demonstrated smaller differences.
Further correlations between DNA methylation levels were calculated separately for each species and cytotype; the results indicate that this adverse relationship between stress and controls was only repeated for diploids of F. amethystina (Fig. 3). However, F. tatrae demonstrated no correlation (0.12) between stress and control conditions in 2018. In addition, tetraploid F. amethystina demonstrated no correlation in methylation level between stress and control in 2018 (− 0.08), nor between stress and controls in 2019 (− 0.09) (Fig. 3).
Factors affecting the level of global DNA methylation. The presence of adverse weather conditions, such as higher insolation and lower air humidity, measured 24 h and not seven days before sampling was asso- www.nature.com/scientificreports/ ciated with an elevated level of DNA methylation (as noted in 2019). In both years, higher temperature and humidity deficiency, indicative of stress conditions, were associated with a relatively higher level of DNA methylation compared to the control periods in the same year ( Supplementary Fig. S1). For both studied species and all cytotypes, an exponential relationship was observed between DNA methylation level and mean insolation (in hours per hour) 24 h before sampling (Fig. 4).
For three of the four analyzed sampling dates, bioclimatic parameters describing the environmental provenance of the individuals were found to be more important for explaining global DNA methylation level than those recorded during the experiment (i.e. phenotypic features of the plants and soil moisture) ( Table 1). The temperature variables were found to give a better explanation of the methylation level than those related to precipitation, and this is particularly clear in the final Model set 3 (Table 1). It is worth noting that one of the most significant variables was the mean diurnal range of temperature (Bio2) in the original locations from where the plants were obtained repeats in models (Table 1).  www.nature.com/scientificreports/ Our results indicate that the level of DNA methylation during stress may be significantly influenced by the climatic provenience of the specimens. The climate parameters from the original location, such as temperature extremes and the range between them, as well as the temperatures characterizing the driest or wettest parts of the year, had the greatest influence. In 2018, Bio2 and mean temperature of driest quarter (Bio9) were most important for explaining DNA methylation level in stress conditions; in 2019, Bio2 was supplemented by a further five variables in stress conditions (Table 1).
Ploidy level was found to be the main factor for explaining methylation level during control conditions in 2018, and Bio 2 for control conditions in 2019. Taken as a whole, the data indicates that the results obtained for the control period in 2018 clearly differ from those of the other periods (Table 1).

Discussion
In our experiment, a higher methylation level was generally observed during stress conditions than in unstressed controls ( Fig. 1), and this is in line with general knowledge on plant epigenetic reactions (cf. 22 ). However, even closely-related species can demonstrate considerable variations in DNA methylation 23 and hence may respond differently to changing environmental conditions 24 . In our case, F. amethystina demonstrated a generally higher  www.nature.com/scientificreports/ mean level of global DNA methylation than its close relative F. tatrae. Although F. amethystina has mixed levels of ploidy, with the tetraploid form having generally higher methylation levels, even the diploid F. amethystina demonstrated higher DNA methylation levels than the diploid F. tatrae, despite having a similar genome size.
However, the present study tested global DNA methylation, without any distinction between CG, CHG and CHH methylation; as such, any conclusions drawn on the differences between the studied species need to be confirmed in further experiments. Furthermore, methylation level can be influenced by numerous other phenomena, such as genome shuffling 25 , genomic shock 26 and subsequent repatterning of expression 27 . It would be easy to attribute the higher level of methylation in tetraploid plants to its larger genome, as in Róis et al. 28 , but research on cytotypes within species indicate that usually no such clear relationship exists (e.g. 29 ). Similar results were also obtained in the present study, with any such changes being dependent on the studied season. Models obtained through step regression (Table 1) indicate that ploidy level plays a key role in explaining the methylation level observed during control conditions in 2018. In 2019, during control conditions, the tetraploid F. amethystina demonstrated lower DNA methylation levels than the diploids (Fig. 1), and the diploids demonstrated weaker vitality than recorded in previous years: most of them were not flowering. These findings suggest that the higher level of DNA methylation observed in diploids may be influenced by the general condition of the plants and their development, and not by ploidy level per se.
It has been reported that DNA methylation patterns may differ between varieties of the same species 30,31 . In the present study, the responses differed between F. amethystina plants and cytotypes, as well as between only diploid F. tatrae (Fig. 2). Similarly, Zheng et al., found varieties of rice to respond to stress in different ways at the DNA methylation level 32 ; these findings were supported by those from further studies on rice indicating that different genotypes, and even tissues, demonstrate differences in cytosine methylation under salinity stress, irrespective of the level of salinity tolerance demonstrated by the genotype 33 .
Our present findings also indicate that, in most cases, plants characterized by a higher methylation level during stress demonstrated smaller differences between stressful and conducive conditions. An exception to this rule was the reaction demonstrated by tetraploids of F. amethystina in 2019, when any significant change in methylation caused by stress was found to return under normal conditions. It could be interpreted that higher methylation levels may cause slower demethylation processes during the transition to non-stress conditions. Wang et al. 13 and Yaish 8 found that in rice subjected to increasing DNA methylation during drought stress, only 70% of the total changes reset to the normal level after returning to non-drought conditions.
The idiosyncratic and variable behavior of the tetraploids of F. amethystina could be explained by the complex structure of the genome, driven in this case by a probable allopolyploid origin 34 . In the case of polyploidy, more suitable DNA methylation modifications can reduce the degree of incompatibility that arises from the presence of two or more genomes in a nucleolus 27 , it also takes part in transposon silencing. However, previous studies of hybrid plants found that while changes in their DNA sequence appeared to be more, or less, additive compared to the parental species, genome methylation and gene expression were not 35 . Another possible explanation for the variable DNA methylation responses displayed by tetraploids of F. amethystina could derive from their wider biogeographical niche and the higher number of original habitats in which they occur; this could potentially drive different adaptations 19 . Hence, further research based on genomic data is required to fully understand the relationship between allopolyploid plasticity and methylation levels.
In plants, DNA methylation seems to be very dynamic and changeable 24 . In our experiment, weather conditions taken 24 h before sampling explain global DNA methylation values better than those taken seven days  www.nature.com/scientificreports/ before sampling. Studies have indicated a change in DNA methylation as a reaction to external factors after seven days of exposure 36 or even after 1 h 37 . In the present study, a particularly strong relationship was observed between DNA methylation level and insolation 24 h before sampling (Fig. 4). This is probably caused by the effect of UV radiation, which is harmful to living organisms, the increase in cell transpiration through the leaves, and the more noticeable drought. Recent research has shown that insolation is a key driver of short-term changes in DNA methylation, and thus the stress experienced by the plant, as well as their physiological changes and adaptation to harmful conditions 38,39 . However, it seems that there is no universal rule governing the DNA methylation response of plants to UV radiation: UV-B radiation elicited DNA demethylation in Artemisia annua 40 and UV-A/B irradiation resulted in minimal changes in DNA methylation in maize 41 . However, the strong and exponential relationship between insolation and global DNA methylation observed in our case is worthy of further studies. Our findings can shed also light on various factors associated with the environmental or biogeographic provenance. Alonso et al. 23 propose a general rule that species with wider geographic ranges tend to demonstrate lower levels of DNA methylation. In our case, the opposite was true: F. tatrae, the narrow endemics, showed a lower level of DNA methylation. However, it is possible that grasses may be the exception to this rule, especially considering the special nature of their epigenome in comparison to eudicots 42 ; indeed, only four of the 279 taxa analyzed by Alonso et al. 23 were grass species.
It has been proposed that populations from different regions and different habitat conditions develop specific methylation patterns; in theory, such patterns help plants optimally match their reaction to the conditions in which they live 22 . Variations in DNA methylation between habitats have been reported in several studies 26,43,44 . Epiloci related to eco-environmental variables, particularly water availability and temperature, have been described for the allotetraploid complex from Dactylorhiza 45 , while DNA methylation differences were reported between vineyards growing in different sub-regions 46 . Our findings suggest that the level of DNA methylation occurring during stress may be significantly influenced by the general climatic provenience of specimens. Furthermore, the greatest influence appeared to be exerted by the general climate parameters in the original location, such as temperature extremes and the range between them, as well as the temperatures characterizing the driest or wettest parts of the year.
In our experiment, the differences observed between plants originating from different regions may be attributed to weather differences between seasons: in 2019, the weather conditions were more stressful than in 2018, and methylation levels appeared to be less dependent on the environmental provenance of individuals. Therefore, the effect of original local adaptation on DNA methylation level in such experiments appears to depend on level and duration of stress. Similar conclusions were stated by Richards et al. 44 , who report that while methylation patterns appear to be partly persistent (induced by original habitat and then maintained), the influence of the bioclimatic parameters of the original locations of the plants are modified by additional elements 44 . It should be noted that most of the regression models calculated in our study poorly describe the level of DNA methylation. It therefore appears that the analyzed variables modify the methylation level rather than independently shaping it.

Materials and methods
Common garden experiment and plant material sampling. Specimens of F. amethystina and F. tatrae were collected during field studies (Supplementary Table S1). These were grown in common garden experimental plots in the Botanical Garden in Lodz, Central Poland (51°45′12.5″N, 19°24′30.2″E). They had been grown together for at least a year before sampling. The experiment was conducted during two growing seasons, i.e. in the years 2018 and 2019. Each year, plant material was taken for total DNA methylation analysis during two periods of different weather conditions: (1) under water and temperature stress and (2) under control (favorable) conditions. The plant material was sampled at the time when plants were fully flowering. In 2018, stress condition samples were taken on 2nd June and controls on 2nd July; in 2019 they were taken on 17th June (stress conditions) and on 12th July (control). Plant material was taken from exactly the same plants during the control and stress conditions, and from precisely the same plants each studied year.
Sampling was performed according to good practice, described by Herrera and Bazaga 47 , with the consideration that variations in methylation level between different organs or developmental stages can differ between plant species 48 . All samples were obtained from the same organ and collected at identical developmental stages: i.e. fully-grown leaves from the middle part of fescue tufts. The sampled leaves were put into zipped plastic bags and frozen (< − 23 °C). Where it was possible, samples were taken from two or three specimens from the same population (biological replicates, Supplementary Table S1). Technical replicates (two per measurement) were performed for all samples during the ELISA test.
Each year, a few days after sampling the material for DNA methylation in control conditions, material was collected for phenotypic measurements.
Plant experiments were performed in accordance with relevant guidelines and regulations; plant material was collected with respective permission if necessary.
Weather conditions during experiment. The day of sampling was chosen according to the ongoing analysis of weather conditions monitored by the meteorological station of the Institute of Meteorology and Water Management (IMGW), located in a similar landscape, 3.1 km away from the experimental plots.
Climatic data was obtained from the Accredited Station of the Institute of Meteorology and Water Management, National Research Institute (pol.: IMGW)-Lodz Lublinek. Daily and hourly data was obtained from the IMGW archives (IMGW website). The following data was used: air humidity (%, accuracy 0.1%), air temperature (°C, accuracy 0.1 °C), humidity deficiency (%, accuracy 0.1%), insolation (hours, accuracy 0.1 h/h). The parameters were calculated for seven days and 24 h prior to collection (Supplementary Fig. S1). Absolute differences in global DNA methylation level between stress and control conditions and variations were calculated in a given year: one calculation for all studied specimens and another for each species and cytotype. The results were analyzed using 'lollipop plots with baseline' in R packages 49 : 'ggplot' and 'Ggally' 50 , and boxplots in Statistica v. 13.3 51 .
The differences in the global DNA methylation level between species, and cytotypes, in a given year, and between sampling conditions (i.e. stress vs. control) were assessed using non-parametric Kruskal-Wallis ANOVA (the data was not normally distributed) in Statistica v. 13.3 51 .
Correlation matrices for the obtained global DNA methylation levels during stress and control conditions in a given year were evaluated separately for each species and cytotype. Pearson's correlation coefficients and their statistical significance were calculated, and correlation matrices were created, using R packages 49 : 'psych' and 'corrplot' .
Factors affecting global DNA methylation. Ploidy level of specimens. Based on earlier data 14,34 , there is no evidence that the F. tatrae were anything other than diploid. For F. amethystina, the ploidy level (diploid or tetraploid) of the individuals growing in the same experimental plots was tested using flow cytometry, in accordance with Rewicz et al. 25 (Supplementary Table S1).
Phenotypic characteristics of specimens. In this study, the following continuous traits were chosen as signs of biological fitness based on previous studies: height of stalks, number of stalks and the number of spikelets on the stalk 52 .
As in the above-mentioned analyses, the material was representative of the species and cytotype: all specimens researched epigenetically were also analyzed morphologically. Comparability of the results was ensured by the fact that the plants were cultivated in the same monitored and controlled environment. From the studied specimens, all stalks were counted. The height and spikelet numbers were measured for five stalks per specimen (Supplementary Tables S2-S5).
Moisture level of the upper horizon of soil during sampling. Measurements of soil moisture level were taken during sampling, directly near the clumps of specimens, up to 5-10 cm below the ground (Supplementary  Tables S2-S5). Values were recorded as % humidity with an accuracy of 0.1% using an SM150 Soil Moisture Kit (delta-t).
Geographic provenance of specimens. Altitude. The altitude of original locations of plants was determined using Google Earth. In addition, a script was developed for assigning the altitude to coordinates of the central points of the population locations (Supplementary Tables S1 and S2-S5).
Bioclimates in original locations. Thirty arc (~ 1 km) resolution raster data was used, incorporating 19 bioclimatic variables from the WorldClim database 53 . Nineteen bioclimatic variables were assigned to each location from the rasters (Supplementary Tables S2-S5). The values were extracted according to the coordinates of locations in the ArcGIS Desktop 9.2: Spatial Analyst tools, Extract values to point tool 54 .

Identification of the most important quantitative factors according to stepwise regression.
Stepwise regression was performed to identify the most important variables. The procedure consisted of iteratively removing predictors in the predictive model, so as to find the subset of variables resulting in the best model. In other words, all predictors (variables from the data set) were initially used in the procedure and then the least contributive predictors were iteratively removed: for the next iteration, the model that yielded the lowest AIC was retained. The procedure stop when removing the next one variable would deteriorate the quality of the model.
For each analyzed period (viz. the stress and control periods for 2018 and 2019-Supplementary Tables S2-S5), the analysis was conducted for the three sets of predictors: Model set (1) soil moisture + phenotypic data + chromosome number; Model set (2) altitude + bioclimatic variables of plant original locations, while Model set (3) included the most important variables from Model set 1 and Model set 2. The limitation of the procedure is that the variables should not be collinear. To identify collinearity among explanatory variables, variance inflation factors (VIF) were used. In our case, a large number of climatic variables were collinear; therefore, we decided not to remove a number of variables from the analysis, but to create data subsets without collinear variables in particular subsets. Following this, a stepwise regression procedure was run in each data subset. In the case of Model set 2, 18 such data subsets were created: each variable appeared in four data subsets, each time in a different variable configuration (18 separate models were used to check whether any of the variables may have an impact on the analyzed phenomenon). The list of variables that were found to be statistically

Data availability
All data generated or analyzed during this study are included in this published article (and its Supplementary Information files). www.nature.com/scientificreports/ Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.