Warm temperature-triggered developmental reprogramming requires VIL1-mediated, genome-wide H3K27me3 accumulation in Arabidopsis Free

Changes in the climate are leading to widespread variations in important environmental factors such as temperature. The average global temperature on Earth has increased by about 0.6°C in the last century (Burnett et al., 2003). Most of this warming has occurred in the past three decades at the rate of 0.15-0.20°C per decade. The increase in global temperature caused by climate change will have a devastating effect on biodiversity, crop yield and, ultimately, on human health (Iizumi et al., 2017; Chen et al., 2018; Tigchelaar et al., 2018; Holbrook et al., 2019; Smale et al., 2019). In response to a high ambient temperature environment, plants undergo thermo-morphogenesis as an adaptive response. Thermo-morphogenesis includes changes in plant architecture, such as elongation of the hypocotyl in seedlings and petiole length, and leaf hyponasty to keep the plant body well above hot ground and to facilitate aeration (Quint et al., 2016). The elevated temperature also affects reproductive development by facilitating floral transition in Arabidopsis (Blázquez et al., 2003). Plants employ sophisticated gene regulatory networks in response to external stimuli, including temperatures, to cope with variable environmental conditions (Lin et al., 2020; Ueda and Seki, 2020). A number of regulatory modules that operate to control the growth and development of plants in response to changing temperatures have been characterized (Lin et al., 2020). Multiple thermo-sensors have also been identified, and it has become apparent that plants employ various routes to reprogram their growth and development in response to temperature fluctuations (Vu et al., 2019; Lin et al., 2020).

Thermo-morphogenesis also includes dramatic changes in transcriptomic profiles in plants (Jung et al., 2016). In Arabidopsis, a histone variant, H2A.Z, is evicted from the histone octamer to be replaced with canonical H2A in response to high ambient temperature at some loci, although the molecular mechanism behind temperature-mediated eviction remains unknown (Kumar and Wigge, 2010; Cortijo et al., 2017; Kumar, 2018; Sakamoto and Kimura, 2018). A component of the SWR1 complex, ARP6, is involved in the deposition of H2A.Z, and loss of ARP6 in Arabidopsis leads to severe developmental and pleiotropic phenotypes (Deal et al., 2005). arp6 mutants exhibit a transcriptomic profile similar to plants exposed to high ambient temperature even at normal ambient temperature (Kumar and Wigge, 2010), including constitutively increased expression of key ambient temperature regulators, such as HEAT SHOCK TRANSCRIPTION FACTOR A2 (HSFA2) and FLOWERING LOCUS T (FT) (Charng et al., 2007; Kumar et al., 2012). A CHD subfamily II ATP-dependent chromatin remodeler, PICKLE (PKL), is another example for which a chromatin regulator controls thermo-morphogenesis (Ho et al., 2013; Zha et al., 2017). PKL has been known to be associated with genes enriched with a repressive histone mark, H3K27me3 (Zhang et al., 2012). pkl mutants display a higher level of H3K27me3 than that of wild type (WT) at the loci of some tested thermo-responsive marker genes, such as INDOLE-3-ACETIC ACID INDUCIBLE 19 (IAA19) and IAA29, and repressed the expression of these genes at high ambient temperature (Zha et al., 2017), resulting in reduced sensitivity in hypocotyl elongation at high ambient temperature.

Vernalization is another example of a temperature response whereby plants become competent to flower at the inductive photoperiod by establishing epigenetic repression of the major floral repressor FLOWERING LOCUS C (FLC) (Bastow et al., 2004; Sung and Amasino, 2004). One important mechanism of epigenetic silencing of FLC is mediated by PRC2. PRC2 becomes enriched at FLC chromatin and increases the level of H3K27me3 by vernalization (De Lucia et al., 2008; Kim and Sung, 2013). The association between PHD finger-containing proteins and PRC2 is important for the efficient deposition of H3K27me3 and gene repression by the Polycomb group (De Lucia et al., 2008). PRC2 biochemically co-purifies with the VERNALIZATION INSENSITIVE3 (VIN3) family proteins, including VIN3, VIL1 (also known as VERNALIZATION 5, VRN5) and VIN3-LIKE 2 (VIL2) (Wood et al., 2006; De Lucia et al., 2008; Kim and Sung, 2013). All VIN3 family proteins, including VIL1, are known to act together with PRC2 to increase repressive histone marks at floral repressor loci in Arabidopsis (Kim and Sung, 2013).

Interestingly, we observed that one of the vernalization (prolonged cold period) insensitive mutants, vil1, displayed hypo-responsive phenotypes to high ambient temperature. VIL1, as a facultative component of PRC2, has been shown to participate in vernalization pathways in Arabidopsis by regulating the expression of FLC (Sung et al., 2006; Greb et al., 2007). VIL1 is also involved in light-dependent repression of growth genes by coordinating PRC2 activity with a photoreceptor, phytochrome B (phyB) (Kim et al., 2021), indicating broader roles of VIL1 in the regulation of developmentally controlled genes.

Here, we show that a regulator of the vernalization (cold temperature) response, VIL1, is also required for the high ambient temperature response. Warm ambient temperature triggers genome-wide accumulation of H3K27me3 in a VIL1-dependent manner and VIL1 is necessary for the transcriptomic change by high ambient temperature. VIL1-mediated H3K27me3 accumulation contributes to the downregulation of genes suppressed by high temperature.

When vil1 mutants were grown under constant temperatures (22°C or 27°C) in short-day (SD) conditions, thermo-induced petiole elongation was attenuated in the vil1 mutants compared with WT at 27°C (Sung et al., 2006) (Fig. 1A,B). In addition, high temperature-mediated acceleration of flowering was significantly attenuated in vil1 mutants (Fig. 1C,D). The expression of FLC was upregulated in vil1 mutants, as previously reported (Kim and Sung, 2013), but was not affected by high temperature either in WT or in vil1 mutants (Fig. S1A). FT, which encodes a key positive regulator of flowering, is induced under high ambient temperature to accelerate flowering (Kumar and Wigge, 2010; Kumar et al., 2012), but FT was not induced in vil1 mutants, consistent with the observed flowering times (Fig. S1B; Fig. 1C,D). A similar hypo-response to high ambient temperature was also observed in vil1 mutant plants grown in long-day (LD) condition (Fig. S1C-G), further demonstrating that VIL1 mediates temperature responses regardless of photoperiod conditions. Our previous study showed that the hypocotyl length of vil1 mutants is longer than that of WT at all temperature conditions. VIL1 synergistically regulates hypocotyl growth with phyB through light-dependent repression of ATHB2 without affecting the expression of a key transcription factor, PHYTOCHROME INTERACTING FACTOR 4 (PIF4), regardless of temperature (Kim et al., 2021). However, when plants were switched from 22°C to 27°C growth conditions, relative hypocotyl lengthening was reduced in vil1 mutants (Fig. S2A). The reduction in hypocotyl elongation response induced by high ambient temperature in vil1 mutants was not because the hypocotyl was already elongated in vil1 mutants at 22°C. Indeed, vil1 mutants were still capable of further elongating their hypocotyl by dark conditions at 22°C (Fig. S2B). Consistent with the hyposensitivity of the vil1 mutant, the transcriptional inducibility of thermo-responsive marker genes, such as YUCCA8 (YUC8) and IAA29, in response to high ambient temperature (Sun et al., 2012; Ma et al., 2016; Hayes et al., 2017; Kim et al., 2020) was also compromised in vil1 mutants (Fig. S2C,D). Two vil1 mutant alleles (vil1-1 and vil1-2) exhibited similar phenotypes and the VIL1 complementation transgenic line (pVIL1:VIL1-myc/vil1-1; Kim et al., 2021) was able to rescue the vil1 mutant phenotypes (Fig. 1; Figs S1, S2), suggesting that VIL1 is required for high ambient temperature responses.

To address the molecular basis of hypo-response of the vil1 mutants towards high ambient temperature, we first performed RNA sequencing (RNA-seq) and defined genes that were differentially regulated in response to constant high ambient temperature (Table S2). We found that 5061 genes were up- or downregulated in response to high temperature in WT, whereas only 910 genes were up- or downregulated in vil1 mutants (Fig. 2A; Table S2). Analysis of differentially expressed genes (DEGs) in WT and vil1 mutants at 27°C versus 22°C revealed that in vil1 mutants the majority of high temperature-suppressed genes were not repressed, and a number of high temperature-induced genes were not activated, as shown in the heatmap and the box plots of high temperature-responsive DEGs (5061) in WT and vil1 (Fig. 2B; Fig. S3A,B).

To identify genes correlated with the hypo-response of vil1 mutants to high temperature, we performed hierarchical clustering analysis and defined 18 clusters (Fig. S3C,D). Furthermore, we performed GO (Gene Ontology) term analysis in all defined clusters (Table S3). We found that of the 18 clusters, cluster 5 and 17 contain genes for which expression was significantly induced by high temperature in WT, but to a lesser extent in vil1 (Fig. S3D). GO term analysis showed that these clusters are enriched with genes related to ‘heat acclimation’, including a crucial transcription factor, HSFA2 (Charng et al., 2007) (Fig. 2C; Table S3). Clusters 10, 12 and 15 were also significantly induced by high temperature in WT, but their expression was already elevated at 22°C in vil1 mutants and further responses to high temperature were compromised compared with those in WT (Fig. S3D). GO term analysis in these clusters enriched terms including ‘regulation of developmental processes’ and ‘response to external biotic stimulus’ (Fig. 2C; Table S3). Overall, clusters 5, 10, 12, 15 and 17 contain VIL1-dependent genes (2723 genes in total) induced by high temperature (Fig. S3E).

By contrast, clusters 7 and 9 contained genes for which expression was strongly suppressed at high temperature in WT, but not in vil1 mutants. Genes in clusters 2, 3 and 11 were also strongly suppressed at high temperature in WT, similar to clusters 7 and 9, but in vil1 they were expressed at lower levels and also showed weakened responses to high temperature (Fig. S3D). Altogether, clusters 2, 3, 7, 9 and 11 contained VIL1-dependent genes (2222 genes in total) suppressed by high temperature (Fig. S3E). Interestingly, GO term analysis of DEGs in this category of genes was strongly enriched with terms related to ‘photosynthesis’ and ‘water deprivation’ (Fig. 2C; Table S3). It has been reported that photosynthetic activity is reduced at high temperature to cope with stress (Djanaguiraman et al., 2018), and the expression of genes induced by water deprivation is also greatly reduced at high temperature to cool down the leaf temperature by increasing transpiration rate (Zandalinas et al., 2016). Collectively, these results suggest that photosynthetic activity and response to water deficit are compromised at high temperature in vil1 mutants. Therefore, our results show that the vil1 mutant is hyposensitive to high ambient temperature in terms of transcriptomic responses as well as of morphological responses.

VIL1 is necessary for the enrichment of H3K27me3 at FLC and its homolog loci during vernalization by forming the PHD-PRC2 complex (Wood et al., 2006; De Lucia et al., 2008; Kim and Sung, 2013). To address whether VIL1 contributes to H3K27me3 accumulation by high temperature, we performed chromatin immunoprecipitation followed by sequencing (ChIP-seq) to detect genome-wide H3K27me3 in WT and vil1 mutants grown either at 22°C or 27°C in SD conditions (Fig. 3). We defined about 8000 peaks and 5000 peaks at both temperatures in WT and vil1, respectively. At 27°C, 4258 loci had differential levels of enriched H3K27me3 compared with those at 22°C in WT (Col-0) (Fig. 3A; Table S4). Remarkably, about 93.5% of changes (3983 out of 4258 loci) were due to an increase in the level of H3K27me3 at 27°C compared with 22°C (Fig. 3B; Table S4), suggesting that a widespread increase in the level of H3K27me3 occurs at many loci in response to high ambient temperature. However, little change in the level of H3K27me3 enrichment in response to high temperature was observed in vil1 mutants (202 loci; Fig. 3A), and only 51 loci out of those 202 loci showed an increase in the level of H3K27me3 (Fig. 3C; Table S4). About 45% (1789 out of 3983) of loci with increased H3K27me3 due to high temperature in WT had lower levels of H3K27me3 in vil1 mutants and the other 25% (979 out of 3983) of loci had accumulated H3K27me3 due to high temperature in a VIL1-dependent manner (Fig. 3D). Therefore, for a total of 70% of the loci high temperature-mediated H3K27me3 accumulation was VIL1 dependent.

Genome-wide accumulation of H3K27me3 was also evident from western blot analysis (Fig. 3E). The level of H3K27me3 was dramatically higher at 27°C compared with 22°C, but there no significant increase was observed in vil1 mutants, confirming that VIL1 mediates the accumulation of H3K27me3 in response to high temperature (Fig. 3E). In addition, the VIL1-dependent accumulation of H3K27me3 at high temperature was observed regardless of photoperiods (Fig. S4A,B). Therefore, our results show that VIL1 contributes, not only to the basal level of H3K27me3, but also to the increased genome-wide accumulation of H3K27me3 in response to high temperature.

To address whether some of the misregulated DEGs in vil1 might have resulted from impaired H3K27me3 accumulation in vil1 at high temperature, we first investigated how many DEGs that were downregulated in response to high temperature in a VIL1-dependent manner were enriched with H3K27me3. Out of all the VIL1-dependent, high temperature-downregulated DEGs (2222), 438 were enriched with H3K27me3. H3K27me3 enrichment was increased in 32% of H3K27me3-enriched, VIL1-dependent, high temperature-downregulated DEGs (140 out of 438 genes; hypergeometric distribution test, P= 7.182797e−17) in response to high temperature (Fig. 3F; Table S4). Interestingly, most of the 140 loci failed to accumulate H3K27me3 in vil1 mutants in response to high temperature (Fig. 3G,H). GO term analysis on these 140 loci revealed that terms related to key developmental regulators, such as ‘regulation of transcription, DNA-templated’ and ‘response to hormone’, were enriched for VIL1-dependent, high temperature-mediated, H3K27me3-accumulated DEGs (Table S4), indicating that H3K27me3 accumulation affects repression by high temperature in this class of genes.

ChIP-seq analysis showed that VIL1 mediates high temperature-induced accumulation of H3K27me3 at high temperature-responsive loci (Fig. 3). Moreover, the level of H3K27me3 increased in a VIL1-dependent manner (Fig. 3E; Fig. S4A,B). To evaluate the molecular relationship between VIL1 and H3K27me3 in high temperature responses, we selected five loci among 140 VIL1-dependent H3K27me3-accumulated DEGs (Fig. 3F), including GLYCINE-RICH PROTEIN 3 (GRP3), COLD-REGULATED 15 B (COR15B), PACLOBUTRAZOL RESISTANCE 6 (PRE6), CELLULOSE SYNTHASE-LIKE G3 (AtCSLG3), which were previously known to be involved in the regulation of plant development by temperature and other stress responses, and a randomly selected gene (AT4G32870) (Fig. 4A). Several studies have demonstrated that GRP3 works as a repressor of cell elongation in cellular stress response and signals (Mangeon et al., 2016). The COR genes respond to low temperature, drought and other abiotic stresses and provide resistance to abiotic stresses (Li et al., 2013; Thalhammer and Hincha, 2014). PRE6 is a transcriptional repressor that negatively regulates the auxin response. PREs have also been reported to regulate the high ambient temperature response by incorporating signals transferred from light, auxin and brassinosteroid pathways (Zheng et al., 2017; Jin and Zhu, 2019). CSLG3 controls plant growth through cell wall synthesis and reorganization by various external stimuli (Yang et al., 2016). The random gene, AT4G32870, encodes a lipid transport superfamily protein.

H3K27me3 ChIP-qPCR and real-time quantitative reverse transcription PCR (qRT-PCR) were performed to verify our H3K27me3 ChIP-seq and RNA-seq results, respectively, for selected loci (Fig. 4B; Fig. S4E,F). For each locus, the P2 region was a VIL1-dependent, H3K27me3-enriched region. Consistent with the ChIP-seq results, H3K27me3 on these five loci accumulated in response to high temperature in a VIL1-dependent manner (Fig. 4B). The expression of these genes was suppressed by high temperature in WT but maintained at high levels in vil1 mutants (Fig. S4E,F). In addition, we evaluated VIL1 enrichment kinetics at target loci in response to high ambient temperature. We confirmed that VIL1 enriches at these loci and, thus, VIL1 directly regulates this set of genes (Fig. 4C). Interestingly, the enrichment of VIL1 was increased by high temperature, especially at the regions where the H3K27me3 was enriched (P2; Fig. 4B), indicating that high temperature triggers the accumulation of VIL1 at target loci to increase the level of H3K27me3. However, the levels of VIL1 transcript and protein were little changed by high temperatures (Fig. S4C,D). Taken together, our results show that VIL1-H3K27me3 accumulation in response to high temperature at these loci is required for proper development of plants in thermo-morphogenesis.

Given that ARP6 affects the chromatin landscape in response to high ambient temperature (Kumar and Wigge, 2010) in an opposite manner to VIL1, we investigated the genetic relationship between VIL1 and ARP6. arp6 mutants phenotypically resemble plants grown under high ambient temperature (Deal et al., 2005; Kumar and Wigge, 2010), including accelerated flowering at normal ambient temperature (22°C) both in LD and in SD conditions (Fig. S5A,B). arp6-1 also displayed longer petioles compared with those of WT in SD conditions at 22°C (Fig. S5C). Interestingly, we found that vil1-1 suppresses the constitutive thermo-morphogenic phenotypes of arp6 mutants in arp6-1 vil1-1 in terms of flowering time and petiole elongation phenotypes. Another allele of vil1 mutant (vil1-2) also suppressed the arp6-1 phenotypes. In addition, the genomic complementing VIL1 construct (VIL1-myc/arp6-1vil1-1) restored arp6 mutant phenotypes in arp6-1 vil1-1 (Fig. S5), confirming the antagonistic effects of mutations in VIL1 on arp6 mutants in ambient temperature responses.

To gain more insight into the antagonistic relationship between VIL1 and ARP6, we additionally performed RNA-seq in seedlings of arp6-1, and arp6-1 vil1-1 grown at 22°C in SD conditions (Table S2) and compared it with the transcriptome of Col-0 or vil1-1 at high ambient temperature by correlation analysis (Fig. S6). Consistent with a previous study (Kumar and Wigge, 2010), the transcriptome of Col-0, WT plants grown at 27°C showed a positive correlation (correlation coefficient, R=0.4) with arp6-1 plants grown at 22°C, indicating that the arp6 mutant at 22°C is similar to WT plants grown at 27°C (Fig. 5A). When we compared the transcriptome of vil1-1 plants grown at 27°C with arp6 plants grown at 22°C, the positive correlation was reduced greatly (R=0.004), consistent with the hypo-response to high temperature of the vil1 mutant (Fig. 5B). Heatmap analysis also supported the results of the hyposensitive response of the vil1 transcriptome to high ambient temperature and the similarity of the arp6 transcriptome to the high-temperature transcriptome of WT (Fig. 5C). Next, we compared transcriptome changes between Col-0, vil1-1, arp6-1 and arp6-1 vil1-1 by examining genes that were up- or downregulated in each genotype in response to high ambient temperature. Consistently, vil1 mutants displayed hyposensitivity to high ambient temperature, and found that arp6 mutants displayed more exaggerated transcriptome changes than those of WT even at normal ambient temperature (Fig. 5D). The exaggerated changes in the transcriptome of arp6 reverted back to the WT levels in arp6 vil1 double mutants (Fig. 5D). In addition, the transcriptome of arp6 vil1 reduced the positive correlation (R=0.4) between the transcriptomes of arp6-1 plants grown at 22°C and of Col-0 plants grown at 27°C to R=0.16 (Fig. 5E). This result is consistent with the reversed phenotype of arp6 vil1 double mutants compared with constitutive thermo-morphogenic arp6 phenotypes. Interestingly, the majority of transcripts that were differentially expressed in arp6 compared with WT reverted back to WT levels in arp6 vil1, resulting in a significant positive correlation (R=0.93) (Fig. 5F), suggesting that VIL1 and ARP6 function antagonistically in transcriptional regulation.

We performed hierarchical clustering of DEGs in Col-0, vil1-1, arp6-1 and arp6-1 vil1-1 to identify genes that are antagonistically affected in arp6 and vil1 mutants (Fig. S7A; Table S4). A large number of genes that belong to either cluster 1 or cluster 2 (each consisting of 3797 and 3840 genes, respectively) showed reversion in arp6 vil1 double mutants (Fig. S7A). GO term analysis identified several categories of genes enriched in clusters 1 and 2 (Table S3). Consistent with the early flowering phenotype of arp6 plants (Fig. S5A,B), cluster 1 enriched GO terms such as ‘positive regulation of flower development’, owing to the presence of RELATIVE OF EARLY FLOWERING 6 (REF6) (Noh et al., 2004), GIGANTEA (GI) (Mizoguchi et al., 2005), CRYPTOCHROME 2 (CRY2) (Guo et al., 1998), LUMINIDEPENDENS (LD) (Lee et al., 1994) and FLOWERING LOCUS K (FLK) (Mockler et al., 2004; Sung et al., 2006; Kim and Sung, 2010) in that cluster. It has been reported that H2A.Z is evicted from the promoter of FT at high ambient temperature, allowing PIF4 binding to FT loci, which promotes FT expression under SD (Kumar et al., 2012). Although FT was not included in any of the defined clusters from our RNA-seq analysis, qRT-PCR analysis showed that increased FT expression in arp6-1 also reverted back to WT levels in arp6-1 vil1-1 (Fig. S1B). Furthermore, cluster 1 enriched a GO term related to ‘positive regulation of growth’ owing to the inclusion of FERONIA (FER) (Guo et al., 2009), MAP KINASE 6 (MPK6) (Meng et al., 2012), ERECTA (ER) (Qu et al., 2017) and LONGIFOLIA 1 (LNG1) (Hwang et al., 2017), suggesting that phenotypic reversion of arp6 in arp6 vil1 may be attributed to transcriptional reversion of these classes of genes in arp6 vil1.

Cluster 2 included genes for which expression was reduced in arp6 but reverted back to WT levels in arp6 vil1 (Fig. S7A). Cluster 2 was enriched with GO terms related to ‘auxin biosynthesis process’ and included eight of the 11 YUCCA (YUC) genes in Arabidopsis (YUC1, YUC3, YUC5, YUC7, YUC8, YUC9, YUC10 and YUC11) (Cao et al., 2019) (Table S4). YUCCA genes encode key enzymes in auxin biosynthesis and have important roles in a wide range of plant developmental processes. The higher order mutants of YUC genes cause severe developmental defects in Arabidopsis (Cao et al., 2019). Low expression of multiple YUC genes in arp6 mutants is correlated with pleiotropic phenotypes of arp6 mutants in various developmental processes (Deal et al., 2005). Cluster 2 also enriched a GO term related to ‘photosynthesis’, and included many genes encoding light harvesting complex subunit proteins (Table S4). This was correlated with the GO term enriched in VIL1-dependent, high temperature-downregulated genes (Fig. S3D,E; 2222 genes) expression of which was suppressed by high temperature in WT, but misregulated in vil1 mutants (Table S3). Indeed, we found that more than 50% of high temperature up- and downregulated, VIL1-dependent DEGs overlapped with clusters 1 and 2, respectively (Fig. S7B,C). Collectively, our result shows that transcriptional changes of a number of developmentally regulated genes in arp6 mutants are reversed in arp6 vil1 double mutants. It should be noted that the number of DEGs (19,033; q<0.05; Table S2) in arp6 mutants compared with WT is far greater than the number of DEGs in vil1 mutants (2599; q<0.05; Table S2) at normal ambient temperature. Therefore, VIL1 appears to specifically function to modulate crucial genes affected by high ambient temperature.

Given the antagonistic relationship between VIL1 and ARP6 in the transcriptome analysis, we tested whether the loss of H2A.Z enrichment can be reversed in vil1 mutants, but most peaks overlapped (Fig. S8A; Table S5), suggesting that there were no dramatic genome-wide changes in H2A.Z distributions in vil1 mutants. We also quantified the relative enrichment of the H2A.Z using deepTools (Ramirez et al., 2014) among those four genotypes, but there was no significant change between WT and vil1, or between arp6 and arp6 vil1 (Fig. S8B). Taken together, these results indicate that overall H2A.Z levels are not affected in vil1 compared with WT at normal ambient temperature.

Next, we examined the changes in H2A.Z level in response to high temperature. We quantified H2A.Z enrichment to identify loci that differentially incorporate H2A.Z in response to high temperature. We identified 1184 loci with different levels of H2A.Z enrichment in WT at high temperature (Fig. 6A). As previously reported (Cortijo et al., 2017), the H2A.Z levels tended to be reduced at 27°C compared with 22°C in most of differentially enriched loci in WT (70%, 829 out of 1184 loci) (Fig. 6B,E). However, differential enrichments of H2A.Z at high temperature at these loci were largely compromised in vil1 mutants (Fig. 6B,F). We identified a similar number of loci (1163 loci) that showed different levels of H2A.Z at high temperature in vil1 mutants, but there was little overlap of these loci with those in WT (152 loci) (Fig. 6C). Furthermore, most loci (∼90%, 1045 out of 1163 loci) that showed different levels of H2A.Z in vil1 mutants accumulated more H2A.Z at high temperature, unlike WT (Fig. 6C). Consistent with the hypo-response to high ambient temperature of vil1 mutants, 80% of the high temperature-mediated H2A.Z-evicted loci in WT (660 out of 829 loci) remained high in vil1 mutants at 27°C without a change in the basal level of H2A.Z at 22°C (Fig. 6D). We also examined the H2A.Z enrichment profile of 660 loci and confirmed that high temperature-mediated H2A.Z eviction in WT at 27°C is significantly impaired in vil1 mutants (Fig. 6G). However, the overall, genome-wide H2A.Z enrichment profiles were not changed in vil1 (Fig. S8B), indicating that the effect of VIL1 on the H2A.Z occupancy is limited to high temperature-responsive H2A.Z eviction.

When we compared 1180 genes for which H2A.Z levels are changed in response to high temperature with high temperature-responsive DEGs (5061, Col-0, 27°C versus 22°C; Fig. 2A), we found that 237 genes were differentially expressed in response to high temperature (Fig. S8C), and enriched with GO terms, including ‘cellular response to heat’ (Table S5), such as HSFA3, and several HEAT SHOCK PROTEIN (HSP) genes. Among these genes, 141 out of 237 DEGs showed a reduction of H2A.Z at high temperature in their loci, and 83% of them (117 out of 141 genes) were VIL1 dependent (Fig. S8D; Table S5). We confirmed that H2A.Z levels at these loci were significantly reduced at 27°C in WT but not in vil1 (Fig. S8E). Interestingly, 117 VIL1-dependent, H2A.Z-evicted DEGs at high temperature displayed a poor gene expression response to high temperature in vil1 mutants, consistent with poor reductions in H2A.Z levels at high temperature in vil1 (Fig. 6H). Eighty out of 117 genes were upregulated at high temperature in a VIL1-dependent manner (Fig. S8F), indicating that VIL1-dependent, H2A.Z-evicted DEGs tend to be upregulated by VIL1 at high temperature in terms of their transcription levels. GO term analysis in 117 VIL1-dependent, H2A.Z-evicted DEGs revealed enrichment with a single term related to ‘cellular process’ (Table S5). This group includes six cell cycle-related genes: MYO-INOSITOL OXIGENASE 1 (MIOX1) (Kanter et al., 2005), CELL DIVISION CYCLE 20.1 (CDC20.1) (Niu et al., 2015), MINICHROMOSOME MAINTENANCE 5 (MCM5) (Shultz et al., 2009), CELLULOSE SYNTHASE A2 (CESA2) (Burn et al., 2002), NOVEL PLANT SNARE 11 (NPSN11) (Zheng et al., 2002) and a R1R2R3 MYB-domain transcription factor, MYB3R-4 (Yang et al., 2021) (Fig. S9A).

We confirmed that high temperature-mediated H2A.Z evictions are compromised in vil1 mutants at several loci by ChIP-qPCR (Fig. S9A-C). We also confirmed that the expression of all cell cycle-related genes in this category are upregulated at high temperature in a VIL1-dependent manner, as shown in our RNA-seq study by qRT-PCR (Fig. S9A; Fig. S10A). However, we could not detect any significant enrichment of VIL1 at these loci by ChIP-qPCR using both VIL1-myc/vil1-1 and VIL1-myc/arp6-1 vil1-1 lines (Fig. S9D), indicating that VIL1 does not directly regulate the expression of these cell cycle-related genes.

There have been several reports that indicate a reciprocal regulation between H3K27me3 and H2A.Z occupancy (Hu et al., 2011; Dai et al., 2017; Carter et al., 2018). However, we found that H3K27me3 was not enriched at these cell cycle-related genes (Fig. S9A). Furthermore, among 658 VIL1-dependent, H2A.Z-evicted loci, 68.2% (449 out of 658) of loci were not enriched with H3K27me3 (Fig. S10B,C), suggesting that most VIL1-dependent H2A.Z eviction occurs independently of H3K27me3. In addition, when we compared VIL1-dependent, H2A.Z-evicted loci with the target loci of CURLY LEAF (CLF) and SWINGER (SWN), two major catalytic subunits of PRC2 (Shu et al., 2019) and with the target loci of FERTILIZATION-INDEPENDENT ENDOSPERM (FIE) (Xiao et al., 2017), a core component of all isoforms of PRC2, and we found that only 5.2% or 7.8% of VIL1-dependent, H2A.Z-evicted genes were also enriched by either CLF/SWN or FIE, respectively (Fig. S10B,C). Taken together, our results suggested that VIL1 indirectly regulates high temperature-mediated H2A.Z eviction at high temperature-responsive loci. Therefore, the genome-wide accumulation of H3K27me3 by VIL1 in response to high temperature occurs in parallel with H2A.Z eviction and the effect of VIL1 on the H2A.Z eviction likely reflects the failure to induce VIL1-dependent, H2A.Z-evicted genes by high temperature.

To explore the antagonism between VIL1 and ARP6 further, we compared 658 VIL1-dependent, H2A.Z-evicted genes with transcriptionally reversed genes in arp6 vil1 (cluster 1 and cluster 2 in Fig. S7A), and found few overlaps, indicating that the antagonistic relationship between VIL1 and ARP6 in transcriptional regulation at normal temperature cannot be explained by the contribution of VIL1 to high temperature-mediated H2A.Z eviction (Fig. S7D). To gain more insights into the role of H3K27me3 in antagonism between VIL1 and ARP6, we performed overlap analysis between transcriptionally reversed genes in arp6 vil1 and H3K27me3-enriched genes (8119) in the Arabidopsis genome. We found that 35% (1343 out of 3840, hypergeometric distribution test, P=9.851008e−58) genes from cluster 2 were enriched with H3K27me3 compared with just 9.2% (351 out of 3797) from cluster 1 (Fig. 6J). Genes in cluster 2 were downregulated in arp6 and reverted in arp6 vil1 (Fig. S7A). We found that H3K27me3 levels of about 50% of H3K27me3-enriched cluster 2 genes (620 out of 1343) were VIL1 dependent (Fig. 6I). This includes six out of the eight YUC genes listed in cluster 2 (Table S4). Therefore, the transcriptional reversion of several key developmental regulators in arp6 vil1 appears to be due to the significant loss of H3K27me3 enrichment in the absence of vil1 and may explain the antagonism between VIL1 and ARP6. It should also be noted that the majority of genes in this class did not show differential H2A.Z deposition in response to high temperature (Fig. S7E). Therefore, VIL1-mediated H3K27me3 deposition at these loci occurs independently of H2A.Z removal at high ambient temperature.

In this study, we report that VIL1 is a regulator in the high ambient temperature response. We found that vil1 mutants are hypo-responsive to high ambient temperature (Fig. 1). Transcriptome analysis revealed that a large fraction of high temperature-responsive genes in WT lost their ability to respond to high ambient temperature in vil1 mutants (Fig. 2). It is intriguing to note that vil1 was originally isolated as a vernalization-insensitive mutant that mediates a low temperature response. Here, we demonstrate that the same protein also plays roles in another part of the temperature spectrum – high ambient temperature.

Polycomb plays important roles, not only in the growth, development and differentiation in eukaryotes (Holec and Berger, 2012), but also in response to abiotic stimuli (Liu et al., 2014; Zong et al., 2022). Mutations in core subunits of PRC2 cause severe developmental defects in plant development (Luo et al., 2000; Shu et al., 2019). Given that VIL1 is a facultative component of PRC2, VIL1 is necessary for the accumulation of H3K27me3 by functioning as a PRC2-associated factor in both low and high temperature responses. Vernalization also triggers the genome-wide accumulation of H3K27me3 in Arabidopsis (Xi et al., 2020). Interestingly, therefore, genome-wide changes in the levels of H3K27me3 induced by both low and high temperatures are biased toward further accumulation of H3K27me3, as high ambient temperature also triggers genome-wide accumulation of H3K27me3 (Fig. 3). In the vernalization response, VIN3, a close homolog of VIL1, is transcriptionally induced by prolonged cold treatment (Sung and Amasino, 2004). The cold-specific induction of VIN3 is a part of a regulatory module to sense and measure the duration of cold temperature by plants (Hepworth et al., 2018). In high temperature responses, various temperature sensors, ranging from phytochromes to RNAs, have been identified (Jung et al., 2016, 2020; Lin et al., 2020). In this study, we showed that VIL1-mediated accumulation of H3K27me3 plays roles in high temperature responses in Arabidopsis. Although H3K27me3 is known as a repressive histone mark, the levels of which are correlated with the silence of genes, global accumulation of H3K27me3 in response to high temperature did not cause global repression of gene expression (Fig. 3F). Given that only about 30% of the direct targets of PRC2 core components (CLF and SWN) are upregulated in the clf;swn double mutants, which show severe developmental defects (Shu et al., 2019), H3K27me3 is likely to create a chromatin landscape that favors gene repression rather than directly driving the transcriptional change. Consistent with that, we observed a correlation between H3K27me3 levels and DEGs in some loci and found that VIL1 accumulation at target loci at high temperature appears to be a part of a regulatory module employed by plants in response to a rise in temperature (Fig. 4).

Our finding that the constitutive thermo-morphogenic phenotype and transcriptome of arp6 mutants is reverted by vil1 mutation in arp6 vil1 mutants (Fig. 5; Fig. S5) implies that VIL1 may modulate the high temperature response by affecting H2A.Z deposition. However, VIL1 does not affect genome-wide distribution of H2A.Z, and H2A.Z levels are completely ARP6 dependent (Fig. S8). Instead, we found that H2A.Z is maintained at high levels in vil1 mutants at high temperature, and is correlated with the misregulation of high temperature-induced genes, including cell cycle-related genes (Fig. 6; Fig. S9). PIF4 mediates high temperature-induced H2A.Z eviction at PIF4 targets together with the INO80 chromatin remodeling complex (Xue et al., 2021), and PIF4 also cooperates with RELATIVE OF EARLY FLOWERING 6 (REF6), which demethylates H3K27me3, to activate high temperature-responsive genes (He et al., 2021). However, VIL1 does not interact with PIF4 (Kim et al., 2021) and indirectly affects high temperature-mediated H2A.Z eviction (Fig. S9). In addition, H2A.Z-evicted loci are predominantly independent of PRC2 and its activity, H3K27me3 deposition (Fig. S10). It has been shown that loci lose H2A.Z during transcriptional activation at high temperature (Sura et al., 2017). Therefore, the failure to evict H2A.Z at VIL1-dependent, H2A.Z-evicted genes may reflect the failure to activate these genes in response to high temperature in vil1 mutants, and that VIL1-mediated H3K27me3 accumulation occurs independently of H2A.Z eviction at high temperature.

In this study, we show that high temperature-mediated VIL1 accumulation at target loci is necessary for H3K27me3 enrichment to downregulate high temperature-suppressed genes and that VIL1 is necessary for the majority of transcriptional changes in response to high temperature. High ambient temperature-mediated accumulation of H3K27me3 by VIL1 at target loci adds another layer to the complicated but sophisticated temperature-sensing machineries in plants. Although our study indicates that H3K27me3 accumulation and H2A.Z eviction at high temperature occur in parallel, it will be interesting to investigate how changes at the level of chromatin is modulated by other known and unknown temperature sensors in plants. In addition, our study expands our understanding of how eukaryotes utilize epigenetic mechanisms under unfavorable environmental conditions. The insights developed from such research will also lead us to develop better strategies to thrive in the era of climate change.

All experiments were carried out in Col-0 ecotype of Arabidopsis thaliana. The Arabidopsis T-DNA lines used in this research were: vil1-1 (SALK_136506), vil1-2 (SALK_140132), arp6-1 (Garlic_599_G03). Plants were grown in LD (16-h light, 8-h dark) and SD (8-h light, 16-h dark) conditions at 22°C or 27°C. pVIL1: VIL1-myc complementation line in vil1-1 mutant background (Kim et al., 2021) was used for VIL1 ChIP-qPCR. Experiment-specific growth conditions are also mentioned in figure legends and/or the Results section. For all experiments, seedlings were germinated and grown on half MS medium, and all samples were harvested at ZT6.

RNA was extracted using TRIzol reagent (Ambion) and 3 µg total RNA was used to synthesize cDNA. After DNase I (Promega) treatment for 30 min at 37°C to remove residual genomic DNA, reverse transcription was performed using M-MLV reverse transcriptase (Invitrogen). The expression levels of target genes were determined by real-time qPCR using specific primer sets (Table S1) and normalized to that of PP2A. Three biological replicates were performed.

For protein analysis, seedling samples were frozen in liquid nitrogen and homogenized well. Total proteins were extracted with urea-denaturing buffer (100 mM NaH₂PO₄, 10 mM Tris-Cl and 8 M urea, pH 8.0, 1 mM PMSF, protease inhibitor cocktails) and the debris was removed by centrifugation (10 min at 12,000 g). The extracted proteins were further denatured by boiling at 100°C for 5 min with 1× SDS sample buffer. Protein levels were detected by immunoblot analysis using an anti-myc antibody (Santa Cruz Biotechnology, c-myc (9E10) X antibody, sc-40 X, 1:10,000 dilution for western blot), anti-H3K27me3 antibody (Millipore, 07-449, 1:500 dilution for western blot), and anti-H3 antibody (Abcam, ab1791, 1:1000 dilution for western blot).

Seedlings were collected and finely ground in liquid nitrogen, then samples were crosslinked in 1% formaldehyde with nuclei isolation buffer and then quenched by adding glycine (final concentration 0.125 M). Sonication was performed by a Bioruptor (30 s on/2 min 30 s off cycles) with high-power output to obtain 200- to 500-bp DNA fragments. Protein–DNA complexes were immunoprecipitated using the following antibodies: anti-myc antibody [Santa Cruz Biotechnology, c-myc (9E10) X antibody, sc-40 X], anti-H3K27me3 antibody (Millipore, 07-449), anti-H3 antibody (Abcam, ab1791) and anti-H2A.Z antibody (Abcam, ab4174). The enrichment of DNA fragments was determined by real-time qPCR using specific primer sets (Table S1). Three biological replicates were performed.

Total RNA was extracted from seedlings (growth conditions described above) using TRIzol (Invitrogen) and treated with DNase I enzyme (Promega) to eliminate traces of genomic DNA. Two biological replicates per each genotype were prepared. Sequencing libraries were prepared with 500 ng total RNA followed by NEBNext Poly(A) mRNA Magnetic Isolation and library preparation using the NEBNext Ultra Directional RNA Library Prep Kit for Illumina (NEB, E7420). The quality of the library was assessed using a Bioanalyzer (Agilent High Sensitivity DNA Assay). All the sequencing was performed on Illumina NextSeq 500 as 35X2 paired end.

For bioinformatics analysis, reads generated by Illumina NextSeq 500 platform (20 million reads) were checked for quality using FastQC. After quality assessment, reads were aligned on the TAIR10 genome using HISAT2. SAM files generated from mapping were then converted into BAM files and sorted using Samtools. Bigwig files for Integrative Genomics Viewer (IGV) visualization were created using deepTools. For gene counting, the bedtools program was used, which generates raw count using BAM files. Raw count was then fed into R for differential gene expression analysis using edgeR and data visualization. For hierarchical clustering analysis, only DEGs were used. Specific analysis details are also provided in figure legends. GO term enrichment analysis was performed using the Gene Ontology website (http://geneontology.org/).

A library was prepared from immunoprecipitated and input samples using NEBNext ChIP-Seq Library Prep Master Mix Set for illumine (NEB, E6200 L). The quality of the library was assessed using Bioanalyzer (Agilent High Sensitivity DNA Assay). Two biological replicates per each genotype were prepared. For bioinformatics analysis, reads generated by Illumina NextSeq 500 platform (15∼25 million reads) were checked for quality using FastQC. After quality assessment, reads were aligned on the TAIR10 genome using Bowtie2. SAM files generated from mapping were then converted into BAM files and sorted using Samtools. Bigwig files for IGV visualization were created using deeptools. MACS2 (Feng et al., 2012) broad peak call was used for peak calling (q<0.01). For differential binding analysis, the DiffBind (Ross-Innes et al., 2012) package was used in R and MA plots were generated using the same package (q<0.05). For H2A.Z enrichment analysis, the deepTools (Ramirez et al., 2014) package was used. First, computeMatrix was generated and then the plotProfile function was used to generate enrichment profiles.

The authors acknowledge the Texas Advanced Computing Center (TACC; http://www.tacc.utexas.edu) at The University of Texas at Austin for providing High Performance Computing resources that have contributed to the research results reported within this paper. We thank Dr Richard Meagher for the anti-H2A.Z antibody and Dr Xiaoyu Zhang for overseeing the H2A.Z ChIP experiment.