O-GlcNAc and phosphorylation modifications on HtL1/FBA10 regulate wheat vernalization for flowering

Domain: biology
Published: 2025-12-27
Editorial outcome: major_revision

Abstract
Vernalization-regulated flowering is vital for wheat yield and geographical distribution, and the diversity of flowering time genes is essential for the breeding of climate-resilient varieties. Sugars have long been recognized in regulating flowering; however, the intrinsic connection between carbohydrate metabolism and vernalization response remains largely unexplored. Here, we identify a fructose 1,6-bisphosphate aldolase (FBA) encoding gene, HtL1/FBA10, as a modulator of heading time variation based on a genome-wide association study utilizing wheat core germplasm collections. Evolutionary analysis shows a decrease in the proportion of haplotype-2 of HtL1, which is linked to delayed flowering, in Chinese and American wheat varieties compared to landraces. Vernalization reduces HtL1/FBA10 phosphorylation levels and increases its O-GlcNAcylation, which in turn enhances its enzymatic activity and facilitates VERNALIZATION1 (VRN1) transcription by regulating histone acetylation at the VRN1 locus. Our findings provide mechanistic insights into the interplay between glucose metabolism and the epigenetic regulation of vernalization in winter wheat.

Introduction
Wheat (Triticum aestivum L.) is the most extensively cultivated staple crop worldwide and constitutes a vital source of caloric intake for the global population1. The growing global population and the effects of climate change present significant challenges to future wheat yield2,3. Wheat is a heterohexaploid species formed through two natural hybridization events during domestication4,5. Flowering timing is vital for grain yield in crops, as it affects environmental conditions during grain filling and impacts reproductive success. Vernalization is crucial for overwintering plants, enabling them to flower after sustained exposure to low temperature6. The gene VRN1 plays a critical role in the genetic network of vernalization response for flowering in wheat, and its expression is contingent upon prolonged cold exposure. VRN1 encodes a MADS-box transcription factor that acts as an integrator to promote flowering7. In temperate cereals like wheat and barley, vernalization-induced expression of VRN1 is intricately regulated by epigenetic modifications8–10. Prior to vernalization, the repressive histone marker H3K27me3 is deposited at the VRN1 locus to inhibit its transcription in barley and wheat. During vernalization, the active markers of H3K4me3 and H3K36me3 gradually accumulate, while H3K27me3 levels decline at VRN18,9. Short-term cold exposure quickly increases histone acetylation at the HvVRN1 locus, activating its transcription through a mechanism similar to the regulation of cold acclimation genes11. However, little is known about the regulatory mechanisms governing histone acetylation at VRN1 in response to prolonged cold for flowering in cereals.

O-linked β-D-N-acetylglucosamine (O-GlcNAc) modification of proteins, reversibly catalyzed by O-GlcNAc transferase (OGT) and O-GlcNAcase (OGA)12, has extensive interplay with phosphorylation13. O-GlcNAcylation and phosphorylation modulate the physiological functions of proteins through synergistic or antagonistic interactions, with the crosstalk between these two post-translational modifications exerting a critical regulatory influence on the pathogenesis of human diseases14–16. O-GlcNAc transferases play a role in the regulation of flowering in plants17,18. In winter wheat, proteins undergo dynamic post-translational phosphorylation and O-GlcNAcylation modifications during the vernalization process19. The vernalization gene VRN-D4, predominantly identified in ancient wheat subspecies from South Asia, represents a duplicated copy of VRN-A1. The RIP-3 region within the first intron of VRN-D4 shows similarity to the binding site of the Arabidopsis GLYCINE-RICH RNA-BINDING PROTEIN 7 (GRP7), which is known to regulate flowering20,21. Three close SNPs in RIP-3 region of VRN-D4 and certain VRN-A1 alleles are linked to vernalization requirements, and variations of these three SNPs influence the binding capacity of GRP2, a wheat homolog of GRP7, to the RIP-3 region20. The binding of GRP2 to RIP-3 is also affected by its post-translational O-GlcNAcylation19, which is catalyzed by TaOGT118. Therefore, O-GlcNAc modification of proteins plays a potential role in vernalization gene-mediated flowering time regulation. Phosphorylation of histone H3 and RNA polymerase II at the FLC locus modulates FLC transcription, affecting flowering transition22,23. Nevertheless, there are a few examples of a single protein exhibiting both O-GlcNAcylation and phosphorylation, referred to as “yin-yang” modifications24, that regulate plant development.

A strong correlation exists between sugar metabolism and epigenetic inheritance in cancer cells. Previous research has demonstrated that metabolites generated through the glycolytic pathway function as substrates or cofactors, playing a significant role in the regulation of epigenetic modification processes25. For instance, acetyl-coenzyme A (acetyl-CoA), serving as the substrate for histone acetylation modification enzymes, influences the levels of histone acetylation mediated by these enzymes within cells26. The dynamic fluctuations in glycolytic activity can impact acetyl-CoA levels, subsequently modulating histone acetylation in both yeast and mammals27,28. Fructose-1,6-bisphosphate aldolase (FBA) facilitates the cleavage of fructose 1,6-bisphosphate (FBP) into dihydroxyacetone phosphate (DHAP) and glyceraldehyde-3-phosphate (G3P), while also catalyzing the reverse reaction to synthesize FBP. FBA plays a crucial role in the glycolytic and gluconeogenic pathways, thereby regulating energy metabolism and carbon flux29. FBA is implicated in various vital biological processes and responds to diverse biotic and abiotic stresses30,31. However, there remains limited understanding of how the glycolysis process integrates FBA activity with specific environmental and developmental signaling in plants.

Here, a genome-wide association study (GWAS) was performed to identify key genetic loci and major genes associated with heading time traits in 91 cultivars of the wheat core germplasm collections. The QTL, quantitative Heading-time Locus 1 (qHtL1), was identified within 639–641 Mb on chromosome 3A, in which a predicted FBA encoding gene, HtL1/FBA10, was identified as a modulator of flowering time. Vernalization leads to reduced phosphorylation and elevated O-GlcNAc modifications of HtL1/FBA10, thereby enhancing its activity. This change subsequently alters histone acetylation levels at the VRN1 locus, promoting flowering in winter wheat.

Results
TraesCS3A02G391100 was identified as the candidate gene of flowering time

Flowering time is a critical factor influencing wheat grain yield and adaptability3. To investigate the genetic variation of heading time among diverse wheat accessions, phenotypic data of heading time of 91 wheat accessions from different regions over three consecutive years (2014, 2015, and 2016), as well as BLUP value heading time over three years, were combined with resequencing data for GWAS analysis32. These accessions comprised 65 modern cultivars, 17 local varieties, and nine elite varieties introduced into China (Supplementary Data 1). The analysis revealed significant association signals across 21 chromosomes, with identified genomic regions containing known flowering-related genes, including FLOWERING LOCUS T 3 (FT3), FRUITFULL-like 3 (FUL3), GIGANTEA (GI), and VRN1 (Supplementary Fig. 1). The major locus qHtL1 associated with flowering was located on chromosome 3A (639–641 Mb) (Fig. 1a, b and Supplementary Fig. 2). A total of 22 genes were located within this region (Supplementary Data 2 and Supplementary Data 3). By analyzing previously published transcriptomic data and publicly accessible databases (http://202.194.139.32/expression/index.html), it was observed that TraesCS3A02G391100, which encodes FBA10, displayed a relatively high expression level during the heading stage (Z-65)4 (Fig. 1c and Supplementary Data 2). Furthermore, TraesCS3A02G391100/FBA10 exhibited the highest expression level in spikes, and TraesCS3A02G391200 and TraesCS3A02G391500 also showed relatively high expression levels in spikes among the 22 genes4 (Supplementary Fig. 3a and Supplementary Data 2). Further analysis of transcriptome data revealed that TraesCS3A02G391100/FBA10 consistently maintained higher expression levels across various developmental stages in spikes33 (Fig. 1d, only genes with FPKM>10 were shown). TraesCS3A02G391200 encodes a protein similar to FLC EXPRESSOR (FLX), which is essential for the establishment of vernalization requirement in Arabidopsis34,35, suggesting its potential role in flowering regulation in wheat. Further haplotype analysis of TraesCS3A02G391200 identified two distinct haplotypes in wheat collections. Phenotypic evaluation revealed no significant differences in heading time among accessions carrying these two haplotypes, indicating the absence of genotype-phenotype association at this locus (Supplementary Fig. 3b). In addition, a haplotype association analysis of TraesCS3A02G391500, which is predicted to encode FBA11, revealed no significant differences in heading time between two haplotypes (Supplementary Fig. 3c). Further analysis revealed two haplotypes for TraesCS3A02G391100/FBA10, with haplotype 2 (Hap-2) being consistent with the Chinese Spring reference genome. Significant differences in heading time were observed between cultivars with these two haplotypes, and cultivars carrying Hap-1 exhibited an earlier heading time (Fig. 1e and Supplementary Fig. 4). Furthermore, the results of real-time fluorescence quantitative PCR (RT-qPCR) showed that the transcript levels of TraesCS3A02G391100/FBA10 in Hap-1 group were generally higher than that in Hap-2 (Supplementary Fig. 5). The SNP (S3A_639082983) in the promoter region of FBA10/HtL1 lies within the predicted binding site for the BBR-BPC transcription factor (Supplementary Data 4), which has been previously validated to regulate gene expression36. Although the amino acid sequences of Hap-1 and Hap-2 remain unchanged, this promoter variation may disrupt BBR-BPC binding, potentially altering FBA10/HtL1 expression. These findings suggested that TraesCS3A02G391100/FBA10/HtL1 (hereafter referred to as HtL1) is the candidate gene potentially associated with this genomic region.

To investigate the geographic distribution and breeding selection of HtL1, the haplotypes were analyzed in 306 worldwide wheat accessions37,38, 355 cultivated varieties bred in different periods39, and 830 Chinese wheat varieties from different breeding periods40. Five main haplotypes (Hap-W1-5) were identified in the 306 worldwide wheat accessions, with Hap-W5 being the predominant haplotype. Hap-W1 predominated in the Middle East and Africa, while Hap-W3 was prevalent in North America (Fig. 1f). Among the 355 landraces and modern cultivars from China and the United States, seven main haplotypes (Hap-CU1-7) were found. Hap-CU2 and Hap-CU7 were the predominant haplotypes. Hap-CU4 and Hap-CU6 were also mainly distributed among American cultivated varieties. Meanwhile, in Chinese cultivated varieties of the 355 accessions, Hap-CU1 and Hap-CU2 showed higher proportions (Fig. 1g). It indicated a priority selection to these haplotypes in China and the United States, possibly related to its genetic adaptability. We identified four haplotypes (Hap-Y1-4) among 830 cultivars released at different breeding stages. Hap-Y4 was the main haplotype, and its proportion gradually decreased from 1961 to 2020 (Fig. 1h).

Using 206 accessions of different ploidy types (comprising 61 diploid accessions, 22 tetraploid accessions and 123 hexaploid accessions), we investigated the variation of HtL1 during wheat domestication. Four main haplotypes (Hap-D1-4) were identified, Hap-D1 was predominantly distributed in wild and domesticated einkorn wheat, while Hap-D2 and Hap-D4 were mainly found in landraces and cultivated varieties, and Hap-D3 was primarily observed in Triticum Urartu (Fig. 1i and Supplementary Fig. 6). The results indicated that different haplotypes may have distinct origin backgrounds. Specifically, Hap-D4 is exclusive to tetraploid and hexaploid wheat, indicating a potential formation during wheat polyploidization. It is hypothesized that Hap-D4 could have originated from Hap-D3 in the process of polyploidization. Previous studies showed that O-GlcNAc modification of proteins regulates sugar metabolism-associated flowering in wheat18,41. Protein modification-omics showed that FBA10/HtL1 undergoes O-GlcNAc and phosphorylation modifications in response to vernalization19. Therefore, HtL1 was selected for further investigation to clarify its role in regulating flowering in wheat.

HtL1 promotes heading in winter wheat

To investigate the biological function of HtL1 in flowering regulation, we constructed overexpression lines of HtL1 by using the constitutive ubiquitin promoter in the winter cultivar KN199 background, thereby enabling the continuous constitutive expression of HtL1-HA in KN199. Meanwhile, we employed the CRISPR-Cas9 gene editing system to knock out all A-, B- and D-homeologs of HtL1 in KN199 background. The target sites for the single guide RNA (sgRNA) were designed within exon 2. In the T0 generation of transgenic plants, only two HtL1 knockout mutants (cr-1 and cr-2) were identified, both of which demonstrated characteristics indicative of fragment deletion (Supplementary Fig. 7). Unfortunately, these two mutants exhibited a lethal phenotype under this condition. A relevant study demonstrates that the knockdown of the gene encoding aldolase A results in lethal outcomes in murine cell lines42. We further generated HtL1 knockdown mutants using RNA interference (RNAi) approach, and a specific 273 bp CDS sequence of HtL1 was selected for precise targeting. Following vernalization, phenotypic analysis of flowering time indicated that the HtL1 overexpression (HtL1-OE) lines exhibited a reduced heading time compared to the KN199 control. Conversely, the HtL1 knockdown mutants displayed a significantly delayed heading time relative to KN199 (Fig. 2a, b and Supplementary Fig. 8). Based on Unique Identifier mRNA Sequencing (UID mRNA-seq) data derived from wheat plumules of KN199 and HtL1-RNAi plants without vernalization treatment, as well as further RT-qPCR analysis of indicated genes at stem elongation stage in HtL1-RNAi and KN199 plants with or without vernalization treatment, we evaluated the impact of HtL1 knockdown on the transcription of other homologous genes, indicating that HtL1 knockdown achieved the anticipated specificity (supplementary Fig. 9). RT-qPCR results revealed that HtL1 expression was upregulated in the HtL1-OE lines but downregulated in the HtL1-RNAi lines compared to KN199 (Fig. 2c). Furthermore, the transcription levels of VRN1 were elevated in the HtL1-OE plants, while decreased in the HtL1-RNAi plants relative to KN199 (Fig. 2d). These findings implied that HtL1 may regulate flowering time in wheat by modulating the expression of VRN1.

Subcellular localization analysis showed that HtL1/FBA10 is located in cytosol (Supplementary Fig. 10). Structural prediction performed using SWISS-MODEL indicated that HtL1/FBA10 forms a tetrameric protein structure, comprising four monomeric units (Fig. 2e, f), which aligns with the structural characteristics of class I aldolases43. To assess the enzyme activity, the HtL1-GST fusion protein was expressed and subsequently purified from Escherichia coli strain BL21 (Supplementary Fig. 11). Enzyme kinetic curves and activity assays collectively demonstrated that HtL1 exhibits fructose-bisphosphate aldolase activity (Fig. 2g, h). Analysis of endogenous aldolase activity in plants demonstrated a significant increase in enzyme activity in HtL1 overexpression lines relative to the KN199, whereas a significant decrease was observed in HtL1 knockdown plants (Fig. 2i). Taken together, these results indicate that the enhancement of aldolase activity of HtL1 facilitate the promotion of flowering in wheat.

TaCDPK13 phosphorylates HtL1 to inhibit its aldolase activity

To explore the regulatory mechanisms governing HtL1 activity, immunoprecipitation-mass spectrometry (IP-MS) was conducted using an anti-HA antibody to identify proteins that potentially interact with HtL1 in cell lysates derived from HtL1-OE plants. The post-translational phosphorylation modification of FBA10/HtL1 was previously observed in response to vernalization19, indicating that the calcium-dependent protein kinase TaCDPK13, among other candidate interaction proteins, may function as a potential interacting partner of HtL1 (Supplementary Data 5). The interaction between HtL1 and TaCDPK13 was validated by a Luciferase Complementation Imaging (LCI) assay in plant. The results showed obvious luminescence signals when Nicotiana benthamiana leaves were co-transformed with HtL1-cLuc and TaCDPK13-nLuc constructs (Fig. 3a). In addition, the interaction between HtL1 and TaCDPK13 is enhanced in the presence of Ca2+ (Fig. 3a). Bimolecular fluorescence complementation (BiFC) analysis revealed robust YFP fluorescence signals in N. benthamiana leaves co-expressing HtL1-cYFP and TaCDPK13-nYFP (Supplementary Fig. 12). Co-immunoprecipitation (Co-IP) further conformed this interaction (Fig. 3b). Collectively, these findings demonstrated that HtL1 physically interacts with TaCDPK13 in vivo and in vitro.

To ascertain whether HtL1 functions as a substrate for TaCDPK13, in vitro experiments were performed. Both the phosphorylation of HtL1 by TaCDPK13 and the autophosphorylation activity of TaCDPK13 were confirmed. The phosphorylation modification levels of HtL1 progressively intensify with the increased incorporation of TaCDPK13 (Fig. 3c), which was further confirmed by mass spectrometry (Supplementary Fig. 13a), and the phosphorylation sites on HtL1 catalyzed in vitro by TaCDPK13 were identified (Supplementary Fig. 13b, c). Furthermore, in vivo constitutively expressed HtL1-HA was immunoprecipitated from HtL1-OE plants using anti-HA antibody for mass spectrometry analysis, and two overlapped phosphorylation sites were identified (Supplementary Fig. 14). These data confirmed that TaCDPK13 is responsible for HtL1 phosphorylation.

To further determine the alterations in phosphorylation levels of HtL1 during vernalization in winter wheat, the phosphorylation status of HtL1-HA from wheat plumules was immunoblotted using anti-phospho Ser/Thr antibody. When maintaining consistent total protein levels across varying durations of vernalization (V0, V7, V14 and V21), the constitutively expressed HtL1-HA protein exhibited accumulation with prolonged cold exposure, whereas its phosphorylation levels showed a significant reduction (Fig. 3d). To assess the impact of TaCDPK13-mediated phosphorylation on the activity of HtL1, HA-tagged HtL1 was expressed in N. benthamiana leaves in the presence or absence of TaCDPK13. The results showed that co-expression of HtL1-HA and GFP-TaCDPK13 led to a reduction in the enzymatic activity of HtL1 when compared to co-expression of HtL1-HA and GFP (Fig. 3e), indicating that the phosphorylation of HtL1 catalyzed by TaCDPK13 inhibits its aldolase activity. The effect of vernalization on in vivo aldolase activity was further evaluated in HtL1-OE wheat plants, revealing that vernalization significantly enhanced aldolase activity (Fig. 3f).

Subsequently, we generated overexpression lines of TaCDPK13 driven by the ubiquitin promoter in KN199 background to investigate the influence of TaCDPK13 on flowering time. Additionally, we identified two mutants, tacdpk13-1 and tacdpk13-2, from the KN9204 EMS mutant library, which induce premature termination of TaCDPK13 protein translation due to base mutations (Supplementary Fig. 15). Under vernalization conditions, TaCDPK13 overexpression lines exhibited a delayed flowering phenotype compared to KN199. In contrast, the two mutants showed an earlier flowering phenotype (Fig. 4a, b). RT-qPCR analysis showed that TaCDPK13 expression was upregulated in the TaCDPK13-OE lines relative to KN199 (Fig. 4c). Furthermore, the expression levels of VRN1 were assessed in the transgenic lines. The transcript abundance of VRN1 was reduced in TaCDPK13-OE plants compared to KN199, whereas it was elevated in tacdpk13 mutants (Fig. 4d, e). These results indicate that TaCDPK13 negatively regulates flowering by inhibiting VRN1 transcription in wheat. To examine the effect of TaCDPK13 on the activity of fructose-1,6-bisphosphate aldolase in wheat, we quantified the overall aldolase activity in TaCDPK13-OE transgenic lines. The results indicated a reduction in activity relative to the KN199 control (Fig. 4f). Incubating cell lysates of TaCDPK13-OE plants with in vitro purified HtL1-GST led to significantly higher HtL1 phosphorylation levels compared to KN199 (Fig. 4g). Together these results suggest that TaCDPK13 negatively regulates flowering by phosphorylating HtL1 to inhibit its activity.

TaOGT1-mediated O-GlcNAc modification of HtL1 enhances its stability

O-GlcNAcylation sites of HtL1 were previously identified through proteomics analysis of protein modifications19, indicating a possible interaction of HtL1 with TaOGT1, an O-GlcNAc transferase reported to regulate flowering in wheat18. First, a yeast two-hybrid assay was performed, which revealed their interaction (Fig. 5a). Furthermore, BiFC, as well as LIC assays, confirmed this interaction in plants (Fig. 5b, c). Previous proteomic studies implied that O-GlcNAcylation of proteins plays a role in regulating vernalization in wheat19. Chemoenzymatic labeling and immunoblot analysis were used to detect the O-GlcNAc modification status of constitutively expressed HtL1 in HtL1-OE lines with or without vernalization. To examine the alterations in O-GlcNAc modification and the overall levels of HtL1 in response to vernalization, we subjected HtL1-OE plants to vernalization treatments with varying durations (V0, V7, V14, and V21). The results demonstrated a progressive accumulation of HtL1, accompanied by an increased proportion of O-GlcNAc modifications of HtL1 (Fig. 5d), following extended vernalization durations.

N-Acetylglucosamine (GlcNAc) can activate intracellular protein O-GlcNAc modification, and O-(2-acetamido-2-deoxy-D-glucopyranosylidene)amino-N-phenylcarbamate (PUGNAc), an inhibitor of O-GlcNAcase (OGA), enhances the levels of intracellular protein O-GlcNAcylation. Alloxan inhibits the activity of OGT, thereby reducing intracellular O-GlcNAc modification levels19,44. During germination, seeds of HtL1-OE plants were treated with combined GlcNAc and PUGNAc, as well as Alloxan alone, and then the constitutively expressed HtL1 in the plumules were detected with anti-HA antibody. Results showed that GlcNAc and PUGNAc together increased HtL1-HA protein levels (Fig. 5e), whereas plants treated with Alloxan showed a reduction in HtL1 protein levels (Fig. 5f). These results indicate that elevated overall O-GlcNAc modification levels resulted in an increase in HtL1 protein levels.

Taken together, these findings suggest that vernalization enhances the O-GlcNAc modification level of HtL1, leading to an elevated abundance of HtL1 by increasing its stability. Additionally, we generated overexpression lines of TaOGT1 driven by the ubiquitin promoter in the KN199 background, and observed that TaOGT1 overexpression plants exhibited early flowering in comparison to KN199 (Supplementary Fig. 16).

HtL1 regulates flowering by mediating histone acetylation modifications at VRN1

Given the role of HtL1 in glycolysis, alterations in HtL1 activity may influence acetyl-CoA production by affecting the final product of the glycolytic pathway, which is related to histone acetylation. Immunoblotting analysis revealed that, compared to KN199, HtL1-RNAi lines had significantly lower levels of H3K14ac and H3K27ac, but no significant change in H3K9ac abundance (Fig. 6a).

Given that HtL1-RNAi lines exhibited a strong delayed flowering phenotype and a significant reduction in mRNA abundance of VRN1, we further explored the correlation between VRN1 transcription and the levels of H3K14ac and H3K27ac at the VRN1 locus. Chromatin immunoprecipitation coupled with quantitative PCR (ChIP-qPCR) was performed to access the enrichment of H3K14ac and H3K27ac in the promoter and first intron regions of VRN1 in A, B and D subgenomes in response to vernalization (Fig. 6b). We concurrently analyzed H3K4me3 levels in these regions. The H3K4me3 modification is enriched in the 5’ promoter and gene body regions of flag leaves and spikelets in winter wheat following prolonged cold exposure, which is consistent with the sustained activation of VRN1 during subsequent warm conditions following vernalization9. We found that under the conditions of a 21-day vernalization treatment, the H3K27ac, H3K14ac, and H3K4me3 levels in I, II, III, IV, V, and VI regions were all decreased in HtL1-RNAi plants compared to KN199 (Fig. 6c–e). To further confirm the relationship between HtL1-catalyzed reaction in the glycolytic pathway and flowering regulation, as well as histone acetylation at VRN1 in response to vernalization, we evaluated the effects of exogenous application of FBP, the substrate for aldolase, on the flowering time of KN199 during the vernalization process. Phenotype observation showed that FBP-treated plants flowered earlier than the untreated control (supplementary Fig. 17a, b). In addition, the H3K14ac and H3K27ac levels at region IV were analyzed, revealing an elevation in both histone modifications (Supplementary Fig. 17c, d), and further RT-qPCR analysis confirmed increased transcription of VRN1 under the same condition (Supplementary Fig. 17e). These findings collectively suggest that HtL1 promotes flowering through elevating H3K14ac and H3K27ac levels at VRN1 to activate its transcription during vernalization process in wheat.

Methods
Plant materials and growth conditions

Wheat (Triticum aestivum L.) cultivar KN199 was used for gene cloning and genetic transformation. Wheat plants were grown at a diurnal temperature regime of 22°C during the day and 18°C at night under long-day conditions (16-h light/8-h dark) in an illuminated plant growth chamber, with a relative humidity of 65%. The seeds were surface-sterilized in 2% NaClO for 20 min, followed by thorough rinsing with flowing water. After imbibition, the seeds were placed on moistened filter papers and allowed to germinate in darkness at 25°C for 3 days, and then were subsequently subjected to a vernalization treatment at 4°C for predetermined durations in complete darkness. Seeds that were germinated and cultivated at 25°C for three days served as the non-vernalized control group (V0). For phenotype observation and tissue harvesting, vernalized or non-vernalized seedlings were transplanted into soil and cultivated in an environmentally controlled greenhouse. Fresh samples were harvested and rapidly frozen in liquid nitrogen and stored at −80°C. Heading time, defined as the date when the primary spike fully emerged from the flag leaf sheath, was statistically analyzed in this study.

GWAS analysis

The GWAS was performed for heading time traits on 91 wheat accessions by GEMMA (Genome-wide Efficient Mixed Model Association) software according to methods previously described37,61, and using the mixed linear model (MLM) analysis.

y=Xα+Sβ+Kμ+e

In this model, y represents phenotype trait; α and β for fixed effects, corresponding to marker effects and non-marker effects respectively; while μ represents unknown random effects. X, S, and K serve as the incidence matrices associated with α, β, and μ, and e is the vector of random residual effects. For constructing the S matrix, the first three principal components (PCs) were utilized through Plink software (version 1.90b6.18). The kinship matrix (K) was established based on the simple matching coefficient matrix. Both the S matrix and the kinship matrix (K) were employed to account for population stratification, thereby correcting for potential population structure effects.

Significant P value thresholds (P<10−5) were set to control the genome-wide type I error rate. Then we enlarged the candidate region to 1000 kb centered on the GWAS signal peak to identity candidate genes. The GWAS results were visualized using the R software package CMplot62, which included the creation of QQ plots and Manhattan plots. Additionally, linkage disequilibrium analysis was performed using the PopLDdecay software63.

Generation of transgenic wheat materials

The full-length coding sequences of HtL1, TaCDPK13, and TaOGT1 were effectively amplified from the winter wheat cultivar KN199. These sequences were then recombined into the pUN1301 vector, which contains the maize ubiquitin promoter to allow constitutive expression of HtL1, TaCDPK13, and TaOGT1, each fused with HA and FLAG tags, respectively.

For the generation of CRISPR/Cas9-mediated FBA10/HtL1 knockout mutants, single guide RNAs (sgRNAs) targeting the second exon were designed using the E-CRISP Design website (http://www.e-crisp.org/E-CRISP.html). Subsequently, the target site sequences were cloned into a CRISPR/Cas9 vector linearized with BsaI (Thermo Fisher Scientific, ER0292) according to the manufacturer’s instructions. RNAi was used to generate HtL1 knockdown mutants. The RNAi vector was constructed to target the second exon of HtL1, a region selected for its minimal sequence identity to other homologous aldolase-encoding genes. All constructs were introduced into Agrobacterium tumefaciens strain EHA105 and then transformed into KN199.

RNA extraction, reverse transcription and RT–qPCR

Total RNA was isolated using the TRIzol reagent (Thermo Fisher Scientific, 15596026), and reverse transcription was performed with the HiScript II Q RT SuperMix by following the manufacturer’s protocol (Vazyme, R223-01). Real-time fluorescence quantitative PCR (RT-qPCR) was conducted on diluted cDNA from three biological replicates using SYBR qPCR Master Mix (Vazyme, Q712-02) in an Mx3000P system. Wheat ACTIN was used to normalize the gene expression levels. The primer sequences for RT-qPCR can be found in Supplementary Data 6.

Library preparation and RNA-Sequencing

Total RNA was isolated using TRIzol Reagent (Thermo Fisher Scientific, 15596026) following standard protocols. Genomic DNA contamination was eliminated through DNase I digestion (extremely England Biolabs, M0303L). RNA purity was verified by measuring the A260/A280 ratio (NanoDropTM OneC, Thermo Fisher Scientific), and RNA integrity was verified via the LabChip GX Touch system (Revvity). RNA concentration was quantified fluorometrically using a Qubit 3.0 system with the Qubit RNA Broad Range Assay kit (Thermo Fisher Scientific, Q10210). mRNA libraries were constructed from 1 μg total RNA using the KC™ Digital mRNA Library Prep Kit (SeqHealth Technology Co., Ltd., Wuhan, China) following the manufacturer’s protocol. To mitigate PCR amplification bias and sequencing errors, the kit incorporates 12-nucleotide unique molecular identifiers (UMIs) for precise labeling of cDNA molecules during first-strand synthesis. Size selection targeting 200–500 bp fragments was performed using magnetic bead-based purification. The enriched libraries were quantified and subjected to PE150 paired-end sequencing on the DNBSEQ-T7 platform (MGI).

Chromatin immunoprecipitation (ChIP)

ChIP assays were conducted following previously established protocols64. Briefly, Total chromatin was fragmented via sonication and subsequently subjected to immunoprecipitation using antibodies against H3K27ac (Abcam, ab4729), H3K14ac (Abcam, Ab52946), or H3K4me3 (Millipore, 07473). The quantification of immunoprecipitated VRN1 chromatin was conducted utilizing qPCR, focusing on the specified regions within the VRN1 locus. The relative abundances of histone modifications were normalized to the input DNA. The primer sequences were provided in Supplementary Data 6.

Immunoblot analysis and Co-IP

Total proteins were extracted from samples and quantified by Bradford assay. The protein samples were boiled at 95°C for 10 min, centrifuged for 10 min, then were separated by 10% SDS-PAGE gels and transferred onto polyvinylidene fluoride membrane (Bio-Rad). The membrane was blocked using 5% (W/V) defatted milk or 5% (W/V) BSA (Bovine Serum Albumin) in 1× Tris-buffered saline containing Tween 20 (TBST, 1× TBS with 0.05% [v/v] Tween-20) at room temperature for 2 h or overnight at 4°C. Subsequently, it was incubated with the primary antibody (anti-Actin [Huaxingbio, HX1843], anti-H3 [Abcam, ab1791], anti-H3K9ac [Abcam, ab218553], anti-CTD110.6 [Cell Signaling Technology, 9875S]) at room temperature for 2 h or at 4°C overnight. Following this, the membrane underwent three washes with 1× TBST for 10 min each time. It was then incubated with an HRP-conjugated secondary antibody (1:5000 dilution) for 1 h at room temperature, followed by three washes with 1× TBST. The chemiluminescence signals were detected using a Tanon-5200 system (Tanon). The band intensities on Western blot images were quantified using ImageJ software.

For Co-IP experiments, the full-length coding sequences of HtL1 and TaCDPK13 were inserted into pCsVMV-HA3-N-1300 and pBI121-GFP vectors, respectively, and then co-transformed into Agrobacterium tumefaciens strain GV3101 and infiltrated into N. benthamiana leaves. After 72 h, total proteins were extracted from the transfected leaves and incubated with anti-GFP magnetic beads (Chromotek, gtma-100) at 4°C for 4 h. Immunoprecipitated proteins were subsequently detected through immunoblot analysis.

Yeast two-hybrid assay

For yeast two-hybrid analysis, the full-length coding sequence of the HtL1 gene was cloned into the pBT3-N vector, while the entire coding region of TaOGT1 was inserted into the pPR3-N vector. To assess protein interactions in yeast, the prey plasmid containing the target gene was co-transformed with the bait construct into the yeast strain MY51 according to the manufacturer’s manual (Clontech, USA). The transformed yeast cells were subsequently cultured on SD medium lacking Trp and Leu (SD/-Trp/-Leu) at 30°C for 2 days. Subsequently, positive transformants were screened on SD/-Trp/-Leu/-His-Ade (20 µg/mL X-gal, 50 mM 3-AT).

Firefly Luciferase Complementation Imaging (LCI) assay

The LCI assay was conducted to detect protein interactions between HtL1 and TaCDPK13 or TaOGT1 in N. benthamiana. Briefly, the full-length coding sequences of HtL1, TaOGT1 and TaCDPK13 were inserted into pCAMBIA1300-cLuc or pCAMBIA1300-nLuc vectors separately. Agrobacterium strain GV3101 containing the specified constructs were co-infiltrated into N. benthamiana leaves. The luminescence images reflecting Luc activity were captured using a CCD imaging system.

BiFC and subcellular localization assays

The full-length coding sequences of HtL1, TaCDPK13 and TaOGT1 were inserted into the pUC-SPYCE and pUC-SPYNE vectors, respectively. Different combinations of constructs were then co-transformed into N. benthamiana leaves. Confocal microscopy (Leica TCS SP5) was employed to observe the fluorescence signals of YFP.

For the HtL1 subcellular localization assay, the coding sequence of HtL1 was cloned into pBI221-GFP and pBI121-GFP vectors, respectively, allowing the expression of HtL1-GFP fusion protein. The control and fusion constructs were independently introduced into Arabidopsis protoplasts by a polyethylene glycol-mediated transient expression system, or infiltrated into N. benthamiana leaves via Agrobacterium (strain GV3101)-mediated transformation. Confocal microscopy (Leica TCS SP5) was employed to observe the fluorescence signals of GFP.

Aldolase activity assay

The full-length coding sequences of HtL1 and TaCDPK13 were individually inserted into pCsVMV-HA3-N-1300 and pBI121-GFP vectors, respectively, and then co-transfected into N. benthamiana leaves. The combinations HtL1-HA and GFP were used as controls. To purify HtL1-HA, total protein was extracted 48 h after transfection using a lysis buffer. The protein sample was incubated with anti-HA magnetic beads (Pierce, 88837) at 4°C overnight, then the HtL1-HA protein was eluted with HA peptide (Pierce, 26184) according to the manufacture’s protocol. The eluent was further purified and concentrated using an Amicon Ultra Centrifugal Filter (Millipore). Aldolase activity of purified HtL1-HA or cell lysates prepared from wheat plumule was determined using a fructose-1,6-bisphosphate aldolase activity assay kit (Solarbio, BC2275) with a microplate reader.

In vitro and in vivo phosphorylation assays

To detect phosphorylation modification of HtL1, the purified recombinant HtL1-GST protein was used as a substrate and incubated with TaCDPK13-GST in a kinase reaction buffer (25 mM Tris-HCl pH 7.5, 5 mM β-Glycerophosphoric acid, 2 mM DL-Dithiothreitol, 0.1 mM Na3VO4, 10 mM MgCl2) at 30°C for 40 min. The reaction was terminated by adding 5× SDS loading buffer and boiling at 95°C for 10 min. A portion of phosphorylated protein was separated by 10% SDS-PAGE and detected using an anti-phospho Ser/Thr antibody65 (1:2000 dilution). The gel stained with Coomassie brilliant blue (CBB) was used as a loading control.

Total protein was extracted from HtL1-OE plants. Whereafter, immunoprecipitation of HtL1 was performed using anti-HA magnetic beads (Pierce, 88837), followed by detection of the samples through 10% SDS-PAGE and immunoblotting with anti-phos Ser/Thr antibody (1:2000 dilution) and anti-HA antibody (Sigma, H9658, 1:5000 dilution).

LC-MS/MS analysis

The fusion proteins HtL1-GST and TaCDPK13-GST, expressed in E. coli strain BL21 (DE3), were purified and incubated in a kinase reaction buffer (25 mM Tris-HCl pH 7.5, 5 mM β-glycerophosphoric acid, 2 mM DL-Dithiothreitol, 0.1 mM Na3VO4, 10 mM MgCl₂). After incubation at 30°C for 40 min, the proteins were subjected to digestion. Specifically, the disulfide bonds were first denatured and reduced by incubation with 5 mM TCEP at 55°C for 30 min. Subsequently, the samples were treated with 25 mM iodoacetamide and incubated in the dark at room temperature for 30 min. Finally, the proteins were digested with approximately 1 μg of sequencing-grade trypsin at 37°C for 16 h, and the resulting samples were analyzed by mass spectrometry. For a detailed experimental procedure, please refer to the supplementary method included in the supplementary information file.

Enzymatic labeling analysis of HtL1 O-GlcNAcylation

Proteins containing O-GlcNAc residues in cell lysates were labeled using enzymatic labeling and biotinylation following the manufacturer’s protocol (Thermo Fisher Scientific, C33368). In brief, total protein was extracted in lysis buffer (20 mM Tris-HCl, pH 7.5, 150 mM NaCl, 10 µM PUGNAc, Protease inhibitor cocktail), and the Click-iT O-GlcNAc Enzymatic Labeling System was utilized to label the cell lysate at 4°C for 24 h. Subsequently, the sample was reacted with biotin alkyne after removing excess UDP-GalNAz, according to the manufacturer’s protocol of the Click-iT Protein Analysis Detection Kit (Thermo Fisher Scientific, C33372). Biotin-labeled samples were then incubated with Streptavidin Magnetic Beads (Beyotime, P2151) overnight at 4°C. Afterwards, anti-HA antibody (Sigma, H9658, 1:5000 dilution) was employed to detect the levels of O-GlcNAc modified HtL1.

Statistics and reproducibility

Data were analyzed using GraphPad Prism software (version 10.4). For comparisons between two groups, an unpaired two-tailed Student’s t-test was applied, while one-way ANOVA was used for comparisons across multiple groups. Significant differences are denoted by asterisks (*P<0.05, **P<0.01, ***P<0.001) or by different lowercase letters (P<0.05). All Western blot experiments were independently repeated three times, with consistent results obtained across replicates. The sample size was not predetermined statistically, and no data were excluded from the analyses. The study was conducted without randomization or blinding of the investigators to group allocation during the experiments and outcome assessment.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "This manuscript studies the vernalization response in wheat. It starts with a genetic analysis and identifies a QTL (qHTL1) associated with heading time after vernalization. They identified the causal gene of qHTL1 as FBA10 that encodes the enzyme fructose 1,6-bisphosphate aldolase. Overexpression of FBA10 caused early flowering whereas knock down caused late flowering, and also either increased or decreased aldolase activity. This correlates with increased or reduced levels of VRN1 mRNA, a promoter of flowering after vernalization. They then argue that post-translational modification of FBA10 is altered during vernalization. Phosphorylation of FBA10 is reduced during vernalization whereas O-GlucNAcylation is increased. These changes are proposed to increase FBA10 protein stability and lead to changes in histone modification at the VRN1 locus that change its expression.\nThis manuscript contains a lot of work going all the way from genetic analysis to detailed biochemical and genomic experiments. However, I was not convinced by several of the conclusions and judged that the data shown were insufficient. I realize the difficulties of working in hexaploid wheat, but perhaps a more focused article that made fewer points but examined them more closely would be more convincing. A smaller point is that the writing needs detailed, perhaps professional, editing in many places.\nDetailed comments:\n1. Resolution of the identification of FBA10. The genetic analysis identifies a large interval including at least 22 genes. The logic for then focusing on FBA10 was suggestive but not totally convincing based on expression level and that the protein is modified during vernalization (lines 122 – 124). Here I wondered what is the difference in sequence between FBA10 alleles, and how does this correlate with gene activity? Also, how many other genes encode FBA10 enzymes in the wheat genome, and how would genetic modifications at this gene affect the total pool of enzyme? Often metabolic enzymes are encoded by genes at multiple loci, making it less clear why modification of a single copy would impact phenotype.\n2. The genetic testing of FBA10 is rather crude (lines 127-129). This is based on very strong overexpression or knock down by RNAi. The knock down could affect other genes as well as the gene of interest and this does not seem to have been tested. It should be possible now even in wheat to make CRISPR alleles of specific genes.\n3. I wondered about the pleiotropy of FBA10 knock down. Since this is an enzyme involved in primary metabolism, do the knock down plants grow more slowly and does this explain at least partly the delay in heading date?\n4. Is the effect on phenotype expected from the effect on gene expression? In Figure 1b-d. I found it surprising that in the OE lines VRN1 expression was increased 4-5x but the heading time phenotype was very small.\n5. I did not find the conclusion that protein modification affected FBA10 protein stability convincing, because I would have liked to see mRNA analysis so that any transcriptional effect could be assessed. Also, I would have liked to see levels of FBA10 protein relative to modified protein assessed in each of the experiments. In Figure 5g, FBA10 protein levels seem to increase during vernalization, but does the level of FBA10 mRNA or is this increase purely due to enhanced stability? In Figure 3d, the level of phosphorylated protein decreases during vernalization, but in this experiment does the level of protein increase as in 5g or not? I would like to see the level of phosphorylated protein relative to the total amount of FBA10 in the same extracts. Same for other modifications.\n6. A small point is that there are many words and phrases requiring editing. In the first two lines of the abstract vernalization is spelt wrongly, and “remains blurry” is not really scientific. In line 32 “involved in two natural hybridization” at least lacks a plural. But there are many examples.", "reviewer_2": "This study describes a candidate gene identified by GWAS on chromosome of 3A from wheat that potentially underlies variation in heading date (the timing of flowering) amongst a collection of predominantly Chinese wheat varieties – qHtL1. The authors propose that a fructose 1,6-diphosphate aldolase 10 is the gene underlying the marker trait associations observed in the region of interest (qHtL1). Then, though biochemical characterisations and analyses of heading date and molecular phenotypes of transgenic plants, the authors provide evidence that biochemical interactions between the candidate gene and a kinase (CDPK13) lead to histone modifications of the VRN1 gene that influence heading date. This is an interesting and worthwhile topic. The authors present a lot of biochemical data to build a mechanistic model to explain at a molecular level the mode of action for qHtL1. Overall the manuscript is clear. There are a number of areas where more information would strengthen the manuscript and support the overall conclusions. The most important point to address is the selection of the candidate gene. At the moment the information used to select the candidate gene was (1) expression data, and (2) statistical association. The weakness here is that all genes in that interval will have a statistical association – the GWAS data tells us that, the genes are all closely linked and will share a common statistical association with the heading date trait. So, the statistical evidence creates circular logic in the current version of the manuscript. Noting also that one other gene in the genetic interval appears to be a known regulator of histone modifications of MADs box genes? (TraesCS3A02G391200 is related to FLX). It would be good to rule out this gene as a candidate. A suggestion is to provide a table summarising sequence expression patterns and natural variation identified in all the candidate genes in that interval. Alternatively, and even more powerful, would be to isolate loss-of-function mutants in the prime candidate gene to provide causal evidence that this gene controls flowering, beyond the transgenic plants used in this study. 1. The abstract could be improved. 2. More information is needed about the source of heading-date data used for GWAS. What were the locations and conditions used for these experiments? 3. It was hard to evaluate the importance of qHtL1 relative to overall variation in heading date. A more comprehensive Manhattan plot should be presented to show the complete picture. 4. Similarly, there needs to be more information provided to support the GWAS analysis. What p-value cut-offs were applied? Was correction of population structure applied? 5. Line 142. More details need to be provided to show how the connection between qHtL1 and CDPK13 was revealed? The method and the data generated are not explained 6. Please check the manuscript for minor errors. For example, in a few instances there are ), at ends of sentences, instead of ). I wish the authors the best of luck with this research.", "reviewer_3": "Yang et al. aim to dissect the mechanisms regulating vernalization for flowering in wheat. Using a genome-wide association study (GWAS), they mapped an important gene associated with flowering time, qHtL1, which encodes fructose 1,6- bisphosphate aldolase (FBA). They showed that overexpression of this gene increased aldolase activity, while RNAi knockdown decreased it, confirming that FBA does indeed have aldolase activity. The overexpression lines showed slightly earlier flowering, while the knockout lines showed delayed flowering, positioning FBA as a key candidate in this region for regulating wheat vernalization using genetic tools. The authors further investigated the regulation of FBA and demonstrated its modulation by CDPK13, a kinase, and TaOGT, which mediates O-GlcNAcylation. Their biochemical data suggest that vernalization decreases the phosphorylation level of FBA, which in turn increases its activity. Overexpression of CDPK13 delayed flowering, while its mutation resulted in an earlier flowering phenotype. In addition, they showed that FBA undergoes O-GlcNAc modification in vivo, with FBA protein levels increasing upon vernalization. Notably, loss of FBA protein affects the state of histone modification on VRN1 chromatin, thereby regulating VRN1 expression. Overall, the paper is interesting. However, several figures have conflicting results and some data need further clarification and quality improvement. In Figure 3d, qHtL1-HA protein levels appear to be consistent across the V0, V7, V14, and V21 stages. However, Figure 5g shows conflicting results, indicating higher qHtL1 protein levels at the V21 stage compared to V0. It is unclear which figure accurately reflects the true levels of qHtL1 and this discrepancy needs to be resolved for consistency. The authors argue that CDPK13 regulates qHtL1, and Figure 3d shows decreased phosphorylation levels using an antibody in a Western blot. It is crucial to identify the exact modification sites in vivo that show this decrease and to confirm whether these sites correspond to those identified in vitro by incubation of TaCDPK13 with qHtL1. If the sites overlap, this would provide stronger evidence that CDPK13 directly modifies qHtL1 and plays a key role in regulating its phosphorylation, thereby affecting qHtL1 activity. A remaining question, however, is how CDPK13 itself is regulated by vernalization. Figure 4g shows an increase in phosphorylation when purified qHtL1-GST is incubated with TaCDPK13-OE4 cell lysate, indicating a slight increase in qHtL1 phosphorylation. However, the quality of the data is insufficient. Western blot quantification is not convincing and requires multiple replicates for validation. In addition, quantification of phosphorylation by alternative methods would increase the reliability of the results. Figure 5 illustrates the regulation of qHtL1 by TaOGT1; however, it is unclear how many replicates were performed to ensure the reproducibility of the data. Furthermore, there seems to be an inconsistency in the interpretation of the data between Figure 5g and Figure 5d. If Figure 5g shows a multi-fold increase in qHtL1-HA protein levels following vernalization treatment, how does this correlate with the consistent input levels shown in Figure 5d? This discrepancy raises questions about the experimental design and normalization methods used. In addition, the number of replicates performed and the reproducibility of these results are not clearly stated, which needs to be addressed to increase confidence in the results. If Figure 5g shows that the qHtL1 level increases several fold between v0 and v21, it is puzzling why the activity of the aldolase level increases only very slightly (900 to 1200?). To prove the model, is it possible to treat the wheat with the product (G3P) of qHtL1 and show that it would affect the histone modification states on VRN1 and affect the flowering time? Additional Points: All figures from Figure 1 through Figure 6 use traditional bar graphs. To improve visual clarity and allow for a better understanding of the actual data distribution, box plots or beeswarm plots should be used."}, "author_response": "Response to Reviewer #1\nReviewer #1 (Remarks to the Author):\nThis manuscript studies the vernalization response in wheat. It starts with a genetic\nanalysis and identifies a QTL (qHTL1) associated with heading time after vernalization.\nThey identified the causal gene of qHTL1 as FBA10 that encodes the enzyme fructose 1,6-\nbisphosphate aldolase. Overexpression of FBA10 caused early flowering whereas knock\ndown caused late flowering, and also either increased or decreased aldolase activity. This\ncorrelates with increased or reduced levels of VRN1 mRNA, a promoter of flowering after\nvernalization. They then argue that post-translational modification of FBA10 is altered\nduring vernalization. Phosphorylation of FBA10 is reduced during vernalization whereas\nO-GlcNAcylation is increased. These changes are proposed to increase FBA10 protein\nstability and lead to changes in histone modification at the VRN1 locus that change its\nexpression.\nThis manuscript contains a lot of work going all the way from genetic analysis to detailed\nbiochemical and genomic experiments. However, I was not convinced by several of the\nconclusions and judged that the data shown were insufficient. I realize the difficulties of\nworking in hexaploid wheat, but perhaps a more focused article that made fewer points but\nexamined them more closely would be more convincing. A smaller point is that the writing\nneeds detailed, perhaps professional, editing in many places.\nDetailed comments:\n1. Resolution of the identification of FBA10. The genetic analysis identifies a large interval\nincluding at least 22 genes. The logic for then focusing on FBA10 was suggestive but not\ntotally convincing based on expression level and that the protein is modified during\nvernalization (lines 122 – 124). Here I wondered what is the difference in sequence\nbetween FBA10 alleles, and how does this correlate with gene activity? Also, how many\nother genes encode FBA10 enzymes in the wheat genome, and how would genetic\nmodifications at this gene affect the total pool of enzyme? Often metabolic enzymes are\nencoded by genes at multiple loci, making it less clear why modification of a single copy\nwould impact phenotype.\nResponse:\nWe appreciate your comment and understand your concerns.\n(1) Resolution of the identification of FBA10.\nIn addition to the expression levels of 22 genes within the significantly associated\ninterval (3A: 639–641 Mb), the association between haplotypes and the phenotype of\nheading date of potential candidate genes was analyzed. Among the 22 genes,\nTraesCS3A02G391200 and TraesCS3A02G391500 exhibited high expression in spike\ntissues, while TraesCS3A02G391100(FBA10) showed significantly higher expression\nlevels across different spike developmental stages compared to other candidate genes in\nthe interval (please see Figure 1a below, Supplementary Figure 3a in the revised version).\nIt is noteworthy that TraesCS3A02G391500 encodes aldolase FBA11.\nthe trait of heading date between two haplotypes of TraesCS3A02G391200\n(Supplementary Figure 3b in the revised version, Figure 1b below). Similarly, the analysis\nfor TraesCS3A02G391500/FBA11 haplotypes also demonstrated no significant differences\nin the heading date trait between the two haplotypes (Supplementary Figure 3c in the\nrevised version, Figure 1c below). In contrast, two haplotypes of FBA10 (Hap1 vs Hap2)\nshowed highly significant differences in heading date phenotypes (Figure 1c in the revised\nmanuscript). On this basis, the posttranslational modification of FBA10 in response to\nvernalization was considered as a clue, indicating a role of FBA10 in environmental and\ngenetic adaptability.\nWe further performed Unique Identifier mRNA Sequencing (UID mRNA-seq) to\nanalyze the transcription levels of TraesCS3A02G391100/FBA10, TraesCS3A02G391200\nand TraesCS3A02G391500/FBA11 in response to vernalization in KN199, and found that\nTraesCS3A02G391100/FBA10 and TraesCS3A02G391200 is upregulated following\nvernalization treatment (V21), whereas TraesCS3A02G391500/FBA11 is downregulated\n(Figure 1d below).\nEvaluating spatiotemporal expression profiles, haplotype-phenotype associations and\nexpression patterns in response to vernalization together, FBA10 was identified as the key\ncandidate gene for subsequent functional investigation.\nFigure 1. Transcription and haplotype analyses of TraesCS3A02G391100/FBA10,\nTraesCS3A02G391200 and TraesCS3A02G391500/FBA11.\n2\n(a) Tissue-specific expression profile of TraesCS3A02G391100/FBA10 and other candidate genes.\nZ: Zadoks decimal code, internationally used to represent developmental stages of cereal crops, with\ndifferent numbers representing distinct developmental stages. The TraesCS3A02G391900 was not\ndetected in any wheat tissues, TraesCS3A02G392500 and TraesCS3A02G392700 exhibited extremely\nlow or no expression across the tested tissues, thus they were not displayed in the figure.\n(b) Heading time analysis between two haplotypes of TraesCS3A02G391200Hap-1 (n=12 accessions) and\nTraesCS3A02G391200Hap-2 (n=58 accessions). The box plots display the interquartile range, comprising\nthe first quartile, median, and third quartile, while the whiskers extend from the minimum to the maximum\nvalues. Two-tailed Student’s t-test for statistical analysis.\n(c) Heading time analysis between two haplotypes of TraesCS3A02G391500Hap-1 (n=13 accessions) and\nTraesCS3A02G391200Hap-2 (n=70 accessions). The box plots display the interquartile range, comprising\nthe first quartile, median, and third quartile, while the whiskers extend from the minimum to the maximum\nvalues. Two-tailed Student’s t-test for statistical analysis.\n(d) The expression patterns of TraesCS3A02G391100, TraesCS3A02G391200 and\nTraesCS3A02G391500 in response to vernalization. KN, wild-type KN199; V0, nonvernalized; V21,\nvernalized for 21days.\n(2) What is the difference in sequence between FBA10 alleles, and how does this\ncorrelate with gene activity?\nResponse：\nWe conducted a detailed comparison of CDS and promoter sequences among the\nthree alleles, TraesCS3A02G391100/FBA10 ， TraesCS3B02G423200/FBA13, and\nTraesCS3D02G383700/FBA15 (see Figure 2 and Figure 3 below). Although these three\ngenes share a higher CDS sequence similarity, their promoter regions display distinct\nvariations in motif composition (Table 1 below). Bioinformatics analysis revealed that the\npromoter region of FBA10 contains significantly more cis-regulatory elements for binding\nof MIKC-type MADS and BBR/BPC transcription factor families, compared to FBA13 and\nFBA15. The MIKC-type MADS-box transcription factor family members (e.g., SOC1, FLC)\ncoordinate signals from the photoperiod pathway and vernalization pathway to regulate\ntranscriptional activation of flowering integrators (e.g., FT, LFY). BBR/BPC transcription\nfactors function in flower development by recruiting Polycomb proteins to mediated\nH3K27me3 at target genes (Ref#1 and Ref#2 below). These analyses suggest that FBA10\nmay function as a central hub in the flowering time regulatory network.\nBased on Unique Identifier mRNA Sequencing (UID mRNA-seq) data, these three\ngenes show distinct transcriptional patterns in response to vernalization. Specifically, the\ntranscription level of FBA10 is the highest both before and after vernalization, with a\nsignificant up-regulation observed following 21 days of vernalization. In contrast,\nvernalization markedly suppresses the transcription of FBA13, whereas FBA15 shows no\nsignificant change (see Figure 4 below), indicating that FBA13 and FBA15 do not\ncontribute to the requirement of high aldose activity for flowering during vernalization\nprocess. Based on the above analysis, the role of FBA10 in promoting flowering under\nvernalization conditions is suggested to be independent of the functions of FBA13 and\nFBA15.\nBinding Studies of BBR/BPC Proteins Lend to Clues in GAGA-Motif Recognition and a Regulatory Role\nin Brassinosteroid Signaling. Frontiers in plant science 10, 466 (2019).\nRef#2: Petrella R, et al. BPC transcription factors and a Polycomb Group protein confine the expression\nof the ovule identity gene SEEDSTICK in Arabidopsis. The Plant journal : for cell and molecular biology\n102, 582-599 (2020).\nFigure 2. CDS sequences alignment of TraesCS3A02G391100/FBA10,\nTraesCS3B02G423200/FBA13 and TraesCS3D02G383700/FBA15.\nTraesCS3B02G423200/FBA13 and TraesCS3D02G383700/FBA15.\nTraesCS3A02G391100/FBA10, TraesCS3B02G423200/FBA13 and TraesCS3D02G383700/FBA15.\nTraesCS3A02G391 TraesCS3B02G423 TraesCS3D02G383\nFamily\n100/FBA10 200/FBA13 700/FBA15\nBBR-BPC 16 0 4\nMIKC-MADS 21 6 12\nC2H2 9 10 5\nNAC 9 5 0\nNin-like 1 0 0\nTrihelix 2 2 0\nDof 7 4 5\nERF 7 51 23\nARF 5 5 1\nB3 2 2 5\nCPP 4 0 0\nG2-like 5 1 1\nMYB 8 3 9\nSBP 2 0 0\nEIL 2 0 1\nHD-ZIP 3 0 7\nC3H 1 1 1\nbZIP 7 1 8\nGATA 1 0 6\nMYB_related 1 0 2\nE2F/DP 0 1 0\nZF-HD 0 1 1\nLBD 0 2 0\nHSF 0 3 1\nTCP 0 1 3\nWRKY 0 0 2\nWOX 0 0 1\nbHLH 0 0 1\nBES1 0 0 4\nTraesCS3B02G423200/FBA13 and TraesCS3D02G383700/FBA15 in response to vernalization.\nNote: Plumules of KN199 with (V21) or without (V0) vernalization were collected for Unique Identifier\nmRNA Sequencing (UID mRNA-seq). V0, Non-vernalized; V21, 21 days of vernalization. Two-tailed\nStudent’s t-test for statistical analysis.\n(3) How many other genes encode FBA10 enzymes in the wheat genome, and how\nwould genetic modifications at this gene affect the total pool of enzyme?\nResponse：\nAs you have noted, the wheat genome contains 21 genes encoding fructose 1,6-\nbisphosphate aldolase. We obtained Unique Identifier mRNA Sequencing (UID mRNA-seq)\ndata from KN199 and FBA10-RNAi plants under V0 conditions, as well as from KN199\nunder both V0 and V21 conditions. The data were used to investigate the expression\nchanges of these 21 alleles in FBA10-RNAi materials and to assess the transcriptional\ndifferences in vernalization response.\nAmong the 21 FBAs identified in wheat, FBA19, FBA20, and FBA21 are classified into\nClass II, while the remaining 18 are distributed across various clusters within Class I, and\nhave distinct subcellular location patterns (see Ref #1 below). Although these FBAs are\ninvolved in distinct biological processes, their collective influence on overall aldolase\nactivity remains largely undetermined. In Escherichia coli, the zinc-dependent class II\nfructose-bisphosphate aldolase (FBA II) accounts for 95 to 100% of the total fructose-\nbisphosphate aldolase activity (Ref #2 below), indicating the different contribution of\naldolases in the total pool of enzyme.\nFBA10 shows higher transcriptional activities across different spike developmental\nstages (Figure 1a above), as well as the whole life cycle (Ref #1 below). In addition, we\nfirstly attempted to generate CRISPR knockout mutants of FBA10. However, the knockout\nmutants displayed lethal phenotypes, indicating a predominant role of FBA10 in regulating\ngrowth and development in wheat. Taken together, it is suggested that FBA10 plays a\ndominant role in the total activity of aldolases pool in wheat.\nthe Fructose-1, 6-bisphosphate Aldolase (FBA) Gene Family in Wheat (Triticum aestivum L.). Front Plant\nSci 8: 1030.\nRef #2：Zgiby, S. M., et al. (2000). Exploring substrate binding and discrimination in fructose1, 6-\nbisphosphate and tagatose 1,6-bisphosphate aldolases. Eur J Biochem 267(6): 1858-1868.\n2. The genetic testing of FBA10 is rather crude (lines 127-129). This is based on very\nstrong overexpression or knock down by RNAi. The knock down could affect other genes\nas well as the gene of interest and this does not seem to have been tested. It should be\npossible now even in wheat to make CRISPR alleles of specific genes.\nResponse:\nThank you for your constructive comments.\nInitially, we engineered CRISPR/Cas9‐based knockout mutants of FBA10. However,\nthe knockout mutants displayed lethal phenotypes and severely impaired developmental\nprocesses, ultimately preventing seed production (also described in the revised\nmanuscript). These findings align with previous literature indicating that FBA knockdown\nleads to lethal effects (Ref #3 and Ref #4 below). In FBA10-RNAi construct, a 273 bp CDS\nsequence exhibiting the highest specificity for FBA10 was selected to minimize off-target\neffects on homologous genes. The effect of FBA10 knockdown on other FBAs encoding\ngenes were evaluated in detail (Supplementary Figure 7 and Figure 6 below).\n1) Kindly refer to the alignment of FBA10 RNAi target sequence and 21 homologous FBA\ngenes sequences corresponding to the target site (Figure 5 below).\n2) Considering that FBA10 serves as a crucial enzyme in basic metabolic processes, its\nknockdown could substantially influence the expression of numerous genes across diverse\npathways. Consequently, we conducted a detailed analysis to assess how FBA10\nknockdown impacts the transcription levels of other homologous FBA genes. Among the\n21 aldolase genes, six genes of FBA11~16 share the greatest sequence homology with\nFBA10. Based on UID mRNA-seq data from KN199 and FBA10-RNAi samples of plumules\nunder V0, the expression levels of FBA11~16 are not significantly affected by FBA10\nknockdown, however, the transcription levels of FBA4, FBA7, FBA8, FBA9, and FBA14\ndecreased in FBA10-RNAi plumules (Figure 6a and b below), in addition to FBA10.\n3) Further analysis showed that vernalization significantly increases the expression of\nFBA10, FBA7, FBA8 and FBA9, decreases FBA4 expression, and has no significant effect\non FBA14 transcription (Figure 6c below), suggesting that FBA4 and FBA14 do not\ncontribute to the enhanced aldolase activity in response to vernalization. Based on the\ndownregulation of FBA7~9 transcription in FBA10-RNAi plumules, we further performed\nqRT-PCR to detect the expression levels of FBA7~9 at the stem elongation stage in\nFBA10-RNAi and KN199 with or without vernalization treatment. We found that none of\nFBA7, FBA8 and FBA9 were downregulated at the stem elongation stage in FBA10-RNAi,\nmoreover, these three genes exhibit similar vernalization response in FBA10-RNAi as that\nin KN199. In contrast, FBA10 keeps a stable significant downregulation under V0 and V21\nconditions (Figure 6c, d and e below).\nTaken together, we suggest that the decreased expression of FBA7, FBA8 and FBA9\n8\nin FBA10-RNAi plumules results from development stage-dependent indirect gene\nexpression regulation, but not off-target effects. These results demonstrate that the 273-\nbp target sequence selected for FBA10 knockdown achieved the anticipated specificity.\nIt is reported that the 21 aldolase in wheat are categorized into three major groups,\nwith FBA10 and FBA7~9 belonging to distinct clusters, exhibiting different subcellular\nlocalization and abiotic stress response patterns. And FBA7~9 are suggested to play roles\nin photosynthesis. Moreover, FBA10 exhibits the highest expression levels during the\ngermination and booting stages among the 21 FBA encoding genes, whereas FBA7~9\nshow the lowest expression levels during these two stages (Ref #1 above).\nIn summary, we propose that FBA10 is the aldolase that plays a pivotal role in\nresponse to vernalization, and the late flowering of FBA10-RNAi transgenic plants under\nvernalization is due to inhibited FBA10 activity. We further suggest that the vernalization-\ninduced changes in the activity and stability of FBA10 are specific, leading to phenotypic\nalterations due to its genetic modification.\nWe hope that our revisions have appropriately addressed your concerns.\nRef #3: Ritterson Lew, C. and D. R. Tolan (2012). Targeting of several glycolytic enzymes using RNA\ninterference reveals aldolase affects cancer cell proliferation through a non-glycolytic mechanism. J Biol\nChem 287(51): 42554-42563.\nRef #4: Zhang, C. S., et al. (2017). Fructose-1,6-bisphosphate and aldolase mediate glucose sensing by\nAMPK. Nature 548(7665): 112-116.\ngenes in wheat.\n(a) Volcano plot showing six significantly downregulated FBA encoding genes in nonvernalized (V0)\nFBA10-RNAi plumules compared to KN199. Nonvernalized wild-type KN199 and FBA10-RNAi plumules\nwere collected for Unique Identifier mRNA Sequencing (UID mRNA-seq).\n(b) Transcriptional levels analysis of FBA4, FBA7, FBA8, FBA9, FBA10, and FBA14 indicates a\ndownregulation, as illustrated in (a). Data are means ± SD, Two-tailed Student’s t-test was used for\nstatistical analysis.\n(c) Expression analysis of FBA4, FBA7, FBA8, FBA9, and FBA14 in response to vernalization in KN199.\nPlumules of KN199 with (V21) or without (V0) vernalization were collected for Unique Identifier\nmRNA Sequencing (UID mRNA-seq). Data are means ± SD, Two-tailed Student’s t-test was used for\nstatistical analysis.\n(d) and (e) qRT-PCR analyses show the relative expression levels of FBA7, FBA8, FBA9, and FBA10\n11\nmRNA during the stem elongation stage in KN199 and FBA10-RNAi plants. KN199 and FBA10-RNAi\nplants were nonvernalized (V0) (d) or vernalized for 21 days (V21) (e), respectively. Leaves were collected\nfor RNA extraction, and the transcript levels of FBAs were normalized to Actin, then normalized to KN199\nplants. Three biological replicates were performed. Data are means ± SD, Two-tailed Student’s t-test was\nused for statistical analysis.\n3. I wondered about the pleiotropy of FBA10 knock down. Since this is an enzyme involved\nin primary metabolism, do the knock down plants grow more slowly and does this explain\nat least partly the delay in heading date?\nResponse：\nThanks for your comment and we understand your concern. To elucidate the effects of\nFBA10 knockdown on wheat growth and development, the morphological phenotypes of\nKN199 and FBA10-RNAi plants at the three-leaf stage was photographed and shown as\nFigure 7a below. Moreover, we analyzed the number of days required for the stem leaves\nto reach the 3rd, 5th, 6th, and 12th growth stages in KN199 and qHtL1-RNAi lines with 21\ndays vernalization. Statistical analysis of three biological replicates revealed no significant\ndifferences in growth rates of stem leaves between KN199 and qHtL1-RNAi plants (Figure\n7b below). These data suggest that the significantly delayed flowering phenotype observed\nin qHtL1-RNAi lines under vernalization conditions is not due to growth retardation from\nFBA10 knockdown.\nFigure 7. Analysis of the influence of qHtL1/FBA10 knockdown on growth rate of KN199.\n(a) Phenotypes of KN199 and qHtL1-RNAi plants at the third-leaf stage. Plants were photographed after\nbeing transplanted to the greenhouse for an additional 9 days of growth, following 21 days of vernalization.\nScale Bar, 10 cm.\n(b) Days required for the indicated stem leaves emerging in KN199 and qHtL1-RNAi wheat plants. The\nstatistical experiment was performed using three biological replicates, with each replicate comprising a\ncount of 12 to 15 plants. Two-tailed Student’s t-test was used for statistical analysis.\n4. Is the effect on phenotype expected from the effect on gene expression? In Figure 1b-\nd. I found it surprising that in the OE lines VRN1 expression was increased 4-5x but the\nheading time phenotype was very small.\nResponse：\nof qHtL1-OE3 plants, demonstrating the significant early flowering phenotype of the qHtL1-\nOE3 line in comparison to KN199 (Figure 8a below).\nVRN1 serves as a key integrator of vernalization and photoperiod signals,\nsubsequently activating the transcription of VRN3/FT. VRN3/FT is transported from the\nleaves to the shoot apical meristem, ultimately facilitating the flowering process in wheat\n(Ref #5 and Ref #6 below). Therefore, we additionally analyzed the expression levels of\nVRN3/FT in leaves of plants after 30 days of growth in the green house, following a\nvernalization treatment of 21 days. The transcriptional level of VRN3/FT in qHtL1-OE3 line\nwas significantly elevated compared to KN199 (Figure 8b below), indicating that qHtL1\noverexpression effectively promotes flowering by activating VRN1, which in turn up-\nregulates VRN3/FT transcription.\nIn addition, we extended vernalization treatment by 30 days for both qHtL1-OE3 and\nKN199 plants, from which we randomly selected fourteen plants and compared their shoot\napical morphology and spike developmental stages after a growth of 28 days in the green\nhouse. We found that, in comparison to KN199, the transition of the stem apex morphology\nto the reproductive growth phase was markedly expedited in qHtL1-OE3 plants compared\nto KN199 (Figure 8c below, and supplementary Figure 6e in the revised manuscript).\nThese findings provide additional evidences that qHtL1 overexpression promotes\nvernalization-induced flowering.\nWe hope these additional experiment results have answered your concern.\nFigure 8. qHtL1 accelerates vernalization-promoted flowering\n(a) Analysis of heading time of qHtL1-OE3 and KN199 plants. Here we show two additional biological\n13\nreplicates, in which KN199 and qHtL1-OE3 plants are grown under a varying soil condition. The box plots\ndisplay the interquartile range, comprising the first quartile, median, and third quartile, while the whiskers\nextend from the minimum to the maximum values, (n=20 plants for each group). Two-tailed Student’s t-\ntest for statistical analysis.\n(b) qRT-PCR analysis shows the relative mRNA abundance of VRN3/FT in KN199 and qHtL1-OE3 plants.\nRNA samples were extracted from leaves of KN199 and qHtL1-OE3 plants at 30 days after transplanted\nin the green house, following 21 days vernalization treatment. The data were normalized to Actin, then\nnormalized to KN199 plants. Independent biological experiments were repeated three times, data are\nmeans ± SD, Two-tailed Student’s t-test for statistical analysis.\n(c) The shoot apical morphology of KN199 and qHtL1-OE3 plants. Plants were vernalized at 4°C for 30\ndays, then transplanted in the green house. After 28 days of growth, a random sample of 14 plants were\nselected for observation and photographic documentation. Scale bars, 3 mm.\nRef #5: Xu, S., & Chong, K. (2018). Remembering winter through vernalisation. Nature Plants, 4(12), 997-\n1009..\nRef #6: Niu, De, et al (2024). \"A molecular mechanism for embryonic resetting of winter memory and\nrestoration of winter annual growth habit in wheat.\" Nature Plants 10.(1): 37-52.\n5. I did not find the conclusion that protein modification affected FBA10 protein stability\nconvincing, because I would have liked to see mRNA analysis so that any transcriptional\neffect could be assessed. Also, I would have liked to see levels of FBA10 protein relative\nto modified protein assessed in each of the experiments. In Figure 5g, FBA10 protein levels\nseem to increase during vernalization, but does the level of FBA10 mRNA or is this\nincrease purely due to enhanced stability? In Figure 3d, the level of phosphorylated protein\ndecreases during vernalization, but in this experiment does the level of protein increase as\nin 5g or not? I would like to see the level of phosphorylated protein relative to the total\namount of FBA10 in the same extracts. Same for other modifications.\nResponse：\nThanks for your comments and concern.\nIn the qHtL1 overexpression construct, qHtL1-HA is driven by the constitutive\nUbiquitin promoter. The presence of the constitutively expressed qHtL1 was verified using specific\nprimers including part of the vector sequence. qRT-PCR analysis revealed no significant\ndifferences in qHtL1 transcription levels across various vernalization periods (V0, V7, V14,\nV21), as illustrated in Figure 9a below, which aligns with expectations. Therefore, the\nobserved variations in qHtL1-HA protein levels during vernalization, as detected by the\nanti-HA antibody, are attributed to differences in protein stability rather than changes in\ntranscription.\nThe original Figure 3d and Figure 5g employed different normalization methods,\nleading to potential ambiguity. The original Figure 3d showed the phosphorylation levels of\nFBA10 across samples with different vernalization durations, with the levels of FBA10-HA\nnormalized for consistency, and the absolute levels of FBA10-HA were not reflected. The\noriginal Figure 5g presented the changes in the absolute FBA10-HA levels over different\n14\ndurations of vernalization treatment, with the levels of total proteins normalized for\nconsistency, showing the gradual increase of constitutively expressed FBA10.\nIn the revised manuscript, the variations in phosphorylation and O-GlcNAc\nmodification levels of constitutively expressed FBA10, relative to the total FBA10-HA at\ndifferent periods of vernalization, were presented utilizing the same normalization method.\nWhen the total protein levels in the V0, V7, V14, and V21 samples were maintained\nconsistently, the total FBA10-HA exhibited a gradual accumulation with prolonged\nvernalization duration, and this trend was consistently observed in Figures 3d and 5d in\nthe revised manuscript. Meanwhile, the phosphorylation levels of FBA10 progressively\ndecreased (Figure 9b below, and Figure 3d in the revised manuscript), while the proportion\nof its O-GlcNAcylation form increased (Figure 9c below, and Figure 5d in the revised\nmanuscript).\nThe dynamic regulation of O-GlcNAc modification of ptoteins is controlled by O-\nGlcNAc transferase (OGT) and O-GlcNAcase (OGA). Upon treatment of cell lysates from\nqHtL1-OE3 plants with PUGNAc, an OGA inhibitor, to elevate the O-GlcNAcylation at total\nprotein levels, the amount of constitutively expressed qHtL1/FBA10 remained relatively\nstable (Figure 5e in the revised manuscript). In contrast, the application of the inhibitor\nAlloxan, which removes O-GlcNAc modifications, resulted in a reduction of FBA10 levels\n(Figure 5f in revised manuscript).\nTaken together, we suggest that O-GlcNAc modification of FBA10 facilitate to enhance\nits stability.\nConsidering the consistent and partly overlapped results, we have replaced the\noriginal Figure 3d with Figure 9b below, and replaced the original Figure 5d with Figure 9c\nbelow. In the revised manuscript, the partially overlapping original Figure 5d and 5g have\nbeen removed.\n(a) The relative mRNA abundance of constitutively expressed qHtL1 in qHtL1-OE3 lines. RNA samples\nwere extracted from plumules of wheat plants vernalized for 0, 7, 14 and 21 days, and transcription levels\nof qHtL1 were analyzed by qRT-PCR. The data were normalized to Actin, then normalized to plants\nvernalized for V0. Independent biological experiments were conducted three times, data are means ± SD,\nTwo-tailed Student’s t-test for statistical analysis.\n(b) Vernalization inhibits phosphorylation of qHtL1. Phosphorylation levels of constitutively overexpressed\nqHtL1 were evaluated during various vernalization periods in qHtL1-OE plants. Total proteins were\nextracted from plumules of wheat plants vernalized at 4°C for 0, 7, 14 and 21 days. Phosphorylation\nmodification of qHtL1 was detected using anti-phosphoSer/Thr antibody, and the overexpressed qHtL1\nwas detected using anti-HA antibody. Actin was probed as a loading control. The phosphorylation levels\nwere firstly normalized to qHtL1, and then normalized to Actin. The relative protein abundance was\nquantified using ImageJ software.\n(c) Vernalization induces accumulation of O-GlcNAcylated qHtL1 in vivo. Chemoenzymatic labelling\napproach was conducted for analyzing O-GlcNAclated qHtL1 in wheat. The qHtL1-OE plants were\nvernalized at 4°C for 0, 7, 14, and 21 days, and then for total protein isolation and chemoenzymatic\nlabelling. Constitutively expressed qHtL1 was detected using anti-HA antibody. The qHtL1 levels were\nnormalized to Actin, and then for quantification of the percentage of O-GlcNAcylated qHtL1 relative its\ntotal amount. The relative protein abundance was quantified using ImageJ software.\n6. A small point is that there are many words and phrases requiring editing. In the first two\n16\nlines of the abstract vernalization is spelt wrongly, and “remains blurry” is not really\nscientific. In line 32 “involved in two natural hybridization” at least lacks a plural. But there\nare many examples.\nResponse：\nThanks for your comments. We have carefully polished our manuscript following the\nsuggestions from you and other reviewers, detailed modifications are highlighted in the\nrevised manuscript with blue fonts, and we hope the revised version is now satisfactory.\nReviewer #2 (Remarks to the Author):\nThis study describes a candidate gene identified by GWAS on chromosome of 3A from\nwheat that potentially underlies variation in heading date (the timing of flowering) amongst\na collection of predominantly Chinese wheat varieties – qHtL1.\nThe authors propose that a fructose 1,6-diphosphate aldolase 10 is the gene underlying\nthe marker trait associations observed in the region of interest (qHtL1).\nThen, though biochemical characterisations and analyses of heading date and molecular phenotypes of\ntransgenic plants, the authors provide evidence that biochemical interactions between the\ncandidate gene and a kinase (CDPK13) lead to histone modifications of the VRN1 gene\nthat influence heading date.\nThis is an interesting and worthwhile topic. The authors present a lot of biochemical data\nto build a mechanistic model to explain at a molecular level the mode of action for qHtL1.\nOverall the manuscript is clear. There are a number of areas where more information would\nstrengthen the manuscript and support the overall conclusions.\nThe most important point to address is the selection of the candidate gene. At the moment\nthe information used to select the candidate gene was (1) expression data, and (2)\nstatistical association. The weakness here is that all genes in that interval will have a\nstatistical association – the GWAS data tells us that, the genes are all closely linked and\nwill share a common statistical association with the heading date trait. So, the statistical\nevidence creates circular logic in the current version of the manuscript.\nNoting also that one other gene in the genetic interval appears to be a known regulator of\nhistone modifications of MADs box genes? (TraesCS3A02G391200 is related to FLX). It\nwould be good to rule out this gene as a candidate.\nResponse：\nThanks for your valuable comments and suggestion.\nThrough genome-wide association study (GWAS) analysis, 22 candidate genes were\nidentified within the chromosomal interval (3A: 639–641 Mb). To investigate their biological\nfunctions in spikelet development, we analyzed transcriptome datasets across multiple\nwheat tissues. The results revealed that TraesCS3A02G391200 and\nTraesCS3A02G391500 exhibited high expression in spike tissues, while\nTraesCS3A02G391100/FBA10 showed significantly higher expression levels across\ndifferent spike developmental stages compared to other candidate genes in the interval\n(Figure 1a below), suggesting its potential involvement in flowering regulation. It is worthy\nto note that TraesCS3A02G391500 encodes aldolase FBA11.\nHaplotype association analysis indicated no significant differences in flowering time\nbetween the two haplotypes of TraesCS3A02G391200 and TraesCS3A02G391500/FBA11,\n18\nrespectively (Figure 1b, c below, Supplementary Figure 3). However, two haplotypes of\nFBA10 (Hap1 vs Hap2) showed highly significant differences in heading date phenotypes\n(Figure 1e in the revised manuscript).\nWe further performed Unique Identifier mRNA Sequencing (UID mRNA-seq) to\nanalyze the transcription levels of TraesCS3A02G391100/FBA10, TraesCS3A02G391200\nand TraesCS3A02G391500/FBA11 in response to vernalization in KN199, and found that\nTraesCS3A02G391100/FBA10 and TraesCS3A02G391200 is upregulated following\nvernalization treatment (V21), whereas TraesCS3A02G391500/FBA11 is downregulated\n(Figure 1d below).\nEvaluating spatiotemporal expression profiles, haplotype-phenotype associations and\nexpression patterns in response to vernalization together, FBA10 was identified as the key\ncandidate gene for subsequent functional investigation.\nFigure 1. Transcription and haplotype analysis of TraesCS3A02G391100/FBA10,\nTraesCS3A02G391200 and TraesCS3A02G391500/FBA11.\n(a) Tissue-specific expression profile of TraesCS3A02G391100/FBA10 and other candidate genes.\nZ: Zadoks decimal code, internationally used to represent developmental stages of cereal crops, with\ndifferent numbers representing distinct developmental stages. The TraesCS3A02G391900 was not\ndetected in any wheat tissues, TraesCS3A02G392500 and TraesCS3A02G392700 exhibited extremely\n19\nlow or no expression across the tested tissues, thus they were not displayed in the figure.\n(b) Heading time analysis between two haplotypes of TraesCS3A02G391200Hap-1 (n=12 accessions) and\nTraesCS3A02G391200Hap-2 (n=58 accessions). The box plots display the interquartile range, comprising\nthe first quartile, median, and third quartile, while the whiskers extend from the minimum to the maximum\nvalues. Two-tailed Student’s t-test for statistical analysis.\n(c) Heading time analysis between two haplotypes of TraesCS3A02G391500Hap-1 (n=13 accessions) and\nTraesCS3A02G391200Hap-2 (n=70 accessions). The box plots display the interquartile range, comprising\nthe first quartile, median, and third quartile, while the whiskers extend from the minimum to the maximum\nvalues. Two-tailed Student’s t-test for statistical analysis.\n(d) The expression patterns of TraesCS3A02G391100, TraesCS3A02G391200 and\nTraesCS3A02G391500 in response to vernalization. KN, wild-type KN199; V0, nonvernalized; V21,\nvernalized for 21days.\nA suggestion is to provide a table summarizing sequence expression patterns and natural\nvariation identified in all the candidate genes in that interval.\nResponse：\nWe appreciate your valuable suggestion. In accordance with your recommendation, we\nhave incorporated Supplemental Data 2, which presents the expression patterns of the 22\ncandidate genes, and Supplemental Data 3, which details the natural variations identified\nwithin all 22 candidate genes in the interval, into the revised manuscript.\nAlternatively, and even more powerful, would be to isolate loss-of-function mutants in the\nprime candidate gene to provide causal evidence that this gene controls flowering, beyond\nthe transgenic plants used in this study.\nResponse：\nWe appreciate your important suggestion and really agree with your comment.\nUnfortunately, the qHtL1 gene mutants that we previously engineered using the\nCRISPR/Cas9 system displayed a lethal phenotype. These findings align with previous\nliterature indicating that FBA knockdown leads to lethal effects (Ref #1 and Ref #2 below).\nFor future functional analyses of other candidate genes, we intend to acquire loss-of-\nfunction mutants, following your recommendations, to further enhance the robustness of\ngene function analysis.\nRef #1: Ritterson Lew, C. and D. R. Tolan (2012). Targeting of several glycolytic enzymes using RNA\ninterference reveals aldolase affects cancer cell proliferation through a non-glycolytic mechanism. J Biol\nChem 287(51): 42554-42563.\nRef #2: Zhang, C. S., et al. (2017). Fructose-1,6-bisphosphate and aldolase mediate glucose sensing by\nAMPK. Nature 548(7665): 112-116.\n1. The abstract could be improved.\nResponse:\nThank you for your valuable suggestion, with which we are in full agreement. Accordingly,\n20\nwe have improved the abstract in the revised version of the manuscript.\n2. More information is needed about the source of heading-date data used for GWAS.\nWhat were the locations and conditions used for these experiments?\nResponse:\nThanks for your suggestions, phenotypic data of heading time of 91 wheat accessions from\ndifferent regions over three consecutive years (2014, 2015, and 2016), as well as BLUP\nvalue heading time over three years, were combined with resequencing data for GWAS.\n3. It was hard to evaluate the importance of qHtL1 relative to overall variation in heading\ndate. A more comprehensive Manhattan plot should be presented to show the complete\npicture.\nResponse:\nThank you for your suggestions. The comprehensive Manhattan plot was presented in\nSupplementary Figure 1 and Figure 2.\n4. Similarly, there needs to be more information provided to support the GWAS analysis.\nWhat p-value cut-offs were applied? Was correction of population structure applied?\nResponse:\nThank you for your suggestions. We have revised this section of the methods pertaining to\nthe GWAS in the revised manuscript. The GWAS was performed for heading time traits on\n91 wheat accessions by GEMMA (Genome-wide Efficient Mixed Model Association)\nsoftware according to methods previously described (Ref #1 and Ref #2 below), and using\nthe mixed linear model (MLM) analysis.\ny = Xα + Sβ + K μ + e\nwhere y represents phenotype; α and β are fixed effects representing marker effects and\nnon-marker effects; and μ represents unknown random effects. X, S, and K are the\nincidence matrices for α, β, and μ, respectively, and e is a vector of random residual effects.\nThe top three PCs were used to build up the S matrix in Plink (version1.90b6.18). The\nmatrix of simple matching coefficients was used to build up the kinship (K) matrix. The S\nmatrix and kinship (K) matrix were performed for population structure correction.\nSignificant P value thresholds (P < 10−5) were set to control the genome-wide type I error\nrate. Then we enlarged the candidate region to 1000 kb centered on the GWAS signal\npeak to identity candidate genes. The GWAS results were visualized using the R software\npackage CMplot (Ref #3 below), which included the creation of QQ plots and Manhattan\nplots. Additionally, linkage disequilibrium analysis was performed using the PopLDdecay\nsoftware. (PopLDdecay: a fast and effective tool for linkage disequilibrium decay analysis\nbased on variant call format files.).\nRef #1: Zhou, X .& Stephens, M (2012). Genome-wide efficient mixed-model analysis for association\nstudies. 482 Nat Genet. 44,821-824(2012).\nRef #2: Liu, Y, et al. (2023). Genetic basis of geographical differentiation and breeding selection for wheat\nplant architecture traits. Genome Biol, 24(1), 114.\nRef #3: Yin, L. et al. (2021). rMVP: A memory-efficient, visualization-enhanced, and parallel-accelerated\n21\ntool for genome-wide association study. Genom Proteom Bioinf. 19, 619-628.\n5. Line 142. More details need to be provided to show how the connection between qHtL1\nand CDPK13 was revealed? The method and the data generated are not explained.\nResponse:\nThanks for your suggestion. In the revised manuscript, we have enhanced the explanation\nregarding the identification of TaCDPK13 as the interacting protein of qHtL1. For further\ndetails, please refer to lines 178-183 in the revised manuscript.\n6. Please check the manuscript for minor errors. For example, in a few instances there\nare), at ends of sentences, instead of).\nI wish the authors the best of luck with this research.\nResponse:\nThanks for your kind comments. We have carefully polished our manuscript following the\nsuggestions from you and other reviewers, detailed modifications are highlighted in the\nrevised manuscript with blue fonts, and we hope the revised version is now satisfactory.\nReviewer #3 (Remarks to the Author):\nYang et al. aim to dissect the mechanisms regulating vernalization for flowering in wheat.\nUsing a genome-wide association study (GWAS), they mapped an important gene\nassociated with flowering time, qHtL1, which encodes fructose 1,6-bisphosphate aldolase\n(FBA). They showed that overexpression of this gene increased aldolase activity, while\nRNAi knockdown decreased it, confirming that FBA does indeed have aldolase activity.\nThe overexpression lines showed slightly earlier flowering, while the knockout lines\nshowed delayed flowering, positioning FBA as a key candidate in this region for\nregulating wheat vernalization using genetic tools.\nThe authors further investigated the regulation of FBA and demonstrated its modulation by\nCDPK13, a kinase, and TaOGT, which mediates O-GlcNAcylation. Their biochemical data\nsuggest that vernalization decreases the phosphorylation level of FBA, which in turn\nincreases its activity. Overexpression of CDPK13 delayed flowering, while its mutation\nresulted in an earlier flowering phenotype. In addition, they showed that FBA undergoes\nO-GlcNAc modification in vivo, with FBA protein levels increasing upon vernalization.\nNotably, loss of FBA protein affects the state of histone modification on VRN1 chromatin,\nthereby regulating VRN1 expression.\nOverall, the paper is interesting. However, several figures have conflicting results and\nsome data need further clarification and quality improvement.\n1. In Figure 3d, qHtL1-HA protein levels appear to be consistent across the V0, V7, V14,\nand V21 stages. However, Figure 5g shows conflicting results, indicating higher qHtL1\nprotein levels at the V21 stage compared to V0. It is unclear which figure accurately reflects\nthe true levels of qHtL1 and this discrepancy needs to be resolved for consistency.\nResponse：\nThank you for your comments and we understand your concerns.\nIn the FBA10-HA overexpression construct, FBA10-HA is driven by the constitutive\nUbiquitin promoter. The anti-HA antibody was employed to monitor changes in the levels\nof constitutively expressed FBA10-HA protein across various vernalization durations.\nThe original Figure 3d and Figure 5g employed different normalization methods,\nleading to potential ambiguity. The original Figure 3d showed the phosphorylation levels of\nFBA10 across samples with different vernalization durations, with the levels of FBA10-HA\nnormalized for consistency, and the absolute levels of FBA10-HA were not reflected. The\noriginal Figure 5g presented the changes in the absolute FBA10-HA levels over different\n23\ndurations of vernalization treatment, with the levels of total proteins normalized for\nconsistency, showing the gradual increase of constitutively expressed FBA10.\nIn the revised manuscript, the variations in phosphorylation and O-GlcNAc\nmodification levels of FBA10, relative to the total FBA10-HA at different periods of\nvernalization, were presented utilizing the same normalization method. When the total\nprotein levels in the V0, V7, V14, and V21 samples were maintained consistently, the total\nFBA10-HA exhibited a gradual accumulation with prolonged vernalization duration.\nMeanwhile, the phosphorylation levels of FBA10 progressively decreased (Figure 1a below,\nand Figure 3d in the revised manuscript), while the proportion of its O-GlcNAcylation form\nincreased (Figure 1b below, and Figure 5d in the revised manuscript). With the progressive\nincrease in vernalization time (V0, V7, V14, V21), there was a significant accumulation of\nthe FBA10-HA protein, and this trend was consistently observed in Figures 3d and 5d in\nthe revised manuscript.\nIn the revised manuscript, we replaced Figure 3d with Figure 1a below，and Figure\n5d with Figure1b below. In the revised manuscript, the partially overlapping original Figure\n5d and 5g have been eliminated.\nFigure 1 Detection of FBA10/qHtL1 phosphorylation and O-GlcNAcylation levels in wheat with\ndifferent vernalization treatments.\n(a) Vernalization inhibits phosphorylation of qHtL1. Phosphorylation levels of constitutively overexpressed\nqHtL1 were evaluated during various vernalization periods in qHtL1-OE plants. Total proteins were\nextracted from plumules of wheat plants vernalized at 4°C for 0, 7, 14 and 21 days. Phosphorylation\nmodification of qHtL1 was detected using anti-phosphoSer/Thr antibody, and the overexpressed qHtL1\nwas detected using anti-HA antibody. Actin was probed as a loading control. The phosphorylation levels\nwere firstly normalized to qHtL1, and then normalized to Actin. The relative protein abundance was\nquantified using ImageJ software.\n(b) Vernalization induces accumulation of O-GlcNAcylated qHtL1 in vivo. Chemoenzymatic labelling\napproach was conducted for analyzing O-GlcNAclated qHtL1 in wheat. The qHtL1-OE plants were\nvernalized at 4°C for 0, 7, 14, and 21 days, and then for total protein isolation and chemoenzymatic\nlabelling. Constitutively expressed qHtL1 was detected using anti-HA antibody. The qHtL1 levels were\nnormalized to Actin, and then for quantification of the percentage of O-GlcNAcylated qHtL1 relative its\ntotal amount. The relative protein abundance was quantified using ImageJ software.\n2. The authors argue that CDPK13 regulates qHtL1, and Figure 3d shows decreased\nphosphorylation levels using an antibody in a Western blot. It is crucial to identify the exact\nmodification sites in vivo that show this decrease and to confirm whether these sites\ncorrespond to those identified in vitro by incubation of TaCDPK13 with qHtL1. If the sites\n24\noverlap, this would provide stronger evidence that CDPK13 directly modifies qHtL1 and\nplays a key role in regulating its phosphorylation, thereby affecting qHtL1 activity. A\nremaining question, however, is how CDPK13 itself is regulated by vernalization.\nResponse：\nThanks for your valuable concerns.\n(1) To further confirm that TaCDPK13 can catalyze the phosphorylation of qHtL1, we\nconducted two supplementary experiments:\n1) Purified qHtL1-GST and TaCDPK13-GST were incubated for an in vitro phosphorylation\nreaction. The reaction products were subsequently analyzed by mass spectrometry (MS)\nto identify the phosphorylation sites on qHtL1. The results indicated that TaCDPK13\ncatalyzed the phosphorylation of qHtL1 at the T33, T35, S330 and S350 sites. Figure 2a\nand b below show the identification of T35 and S350 by MS analysis (also shown as\nsupplementary Figure 11 in the revised version).\n2) We used anti-HA antibody to enrich qHtL1-HA from the cell lysate of nonvernalized\nqHtL1-OE3 plants, and then performed mass spectrometry analysis to detect\nphosphorylation sites on qHtL1, confirming that T35 and S350 sites undergo\nphosphorylation in vivo (Figure 3a, b below and supplementary Figure 12 in the revised\nversion).\nBoth methods above identified overlapped phosphorylation sites of T35 and S350 on\nqHtL1. Thus, our findings substantiate that TaCDPK13 directly catalyze the\nphosphorylation of qHtL1 in vitro and in vivo.\nFurthermore, phosphorylation of qHtL1 at site S350 was also identified in our\npreviously published article with a study on modification proteomics for vernalization\nresponsiveness in winter wheat (Ref #1 below).\nRef1 #: Xu, S, et al. (2019). The protein modifications of O-GlcNAcylation and phosphorylation mediate\nvernalization response for flowering in winter wheat. Plant physiology, 180(3), 1436-1449\nat T35 and S350 sites of qHtL1 expressed in E. coli.\n(a) Identification of phosphorylation at T35 site on qHtL1 using LC-MS. Peptide sequence including T35\nsite, GILAADESTGT(35)IGKR.\n(b) Identification of phosphorylation at S350 site on qHtL1 using LC-MS. Peptide sequence including S350\nsite, ANSEATLGTYKGDATLGEGASES(350)LHVK.\nNote: Fusion proteins of qHtL1-GST and TaCDPK13-GST were expressed and purified in vitro.\nSubsequently, they underwent an in vitro phosphorylation reaction, followed by analysis using liquid\nchromatography-mass spectrometry (LC-MS).\nat T35 and S350 sites on qHtL1 in wheat plants.\n(a) Identification of phosphorylation at T35 site on qHtL1 using LC-MS. Peptide sequence including T35\nsite, GILAADESTGT(35)IGK.\n(b) Identification of phosphorylation at S350 site on qHtL1 using LC-MS. Peptide sequence including S350\nsite, ANSEATLGTYKGDATLGEGASES(350)LHVK.\nTotal proteins were extracted from nonvernalized qHtL1-OE3 wheat plants. qHtL1-HA was\nimmunoprecipitated using an anti-HA antibody, and its phosphorylation modification sites were analyzed\nby LC-MS.\n(2) How CDPK13 itself is regulated by vernalization?\nWe analyzed the transcription levels of TaCDPK13 during different vernalization periods.\nThe data showed that the transcription levels of TaCDPK13 progressively increased with\nthe prolonged vernalization duration (Figure 4 below). As the phosphorylation modification\nlevel of qHtL1 progressively declines during vernalization, we hypothesize that TaCDPK13\nmay also be involved in other signaling pathways. This hypothesis requires further\n27\nexperimental validation, such as examining the binding affinity between TaCDPK13 and\nqHtL1 during different vernalization periods.\nFigure 4. The relative mRNA abundance of TaCDPK13 at different vernalization\ndurations in KN199.\nRNA samples were extracted from plumules of wheat plants vernalized for 0, 7, 14 and 21 days, and\ntranscription levels of TaCDPK13 were analyzed by qRT-PCR. The data were normalized to Actin, then\nnormalized to plants vernalized for V0. Independent biological experiments were conducted three times,\ndata are means ± SD, Two-tailed Student’s t-test for statistical analysis.\n3. Figure 4g shows an increase in phosphorylation when purified qHtL1-GST is incubated\nwith TaCDPK13-OE4 cell lysate, indicating a slight increase in qHtL1 phosphorylation.\nHowever, the quality of the data is insufficient. Western blot quantification is not convincing\nand requires multiple replicates for validation. In addition, quantification of phosphorylation\nby alternative methods would increase the reliability of the results.\nResponse：\nThank you for your comments. To address your concern, we conducted two additional\nbiological replicates, which demonstrated an increase in qHtL1 phosphorylation level and\nenhanced data quality, and the representative data was shown as Figure 5 below. We\nreplaced the original Figure 4g with Figure 5 below with a higher quality in the revised\nmanuscript.\nTo further confirm that TaCDPK13 can catalyze the phosphorylation of qHtL1, we\nconducted a supplementary experiment to identify phosphorylation sites and quantify\nphosphorylation levels of qHtL1. Purified qHtL1-GST and TaCDPK13-GST were incubated\nfor an in vitro phosphorylation reaction. The reaction products were subsequently analyzed\nby mass spectrometry to identify the phosphorylation sites on qHtL1. The data indicated\nthat TaCDPK13 catalyzes the phosphorylation of qHtL1 at the S330, S348, S350, Y18, and\nT35 sites. Furthermore, MS results demonstrated that increasing the amount of TaCDPK13\nin the reaction system correspondingly enhances the degree of phosphorylation at these\nsites (Figure 6 below, and Supplementary Figure 11a).\nWestern blot together with LC-MS analyses demonstrated that TaCDPK13 is\n28\nresponsible for qHtL1phosphorylation.\nFigure 5. Overexpressed TaCDPK13 in wheat plants phosphorylates purified qHtL1.\nIn vivo overexpressed TaCDPK13 increases the phosphorylation level of in vitro purified qHtL1. Total\nproteins were isolated from plumules of TaCDPK13-OE lines and KN199 without vernalization, followed\nby incubation with purified qHtL1-GST for subsequent phosphorylation analysis. Actin was used as a\nloading control. The relative protein abundance was quantified using ImageJ software.\nFigure 6. Detection of variation of qHtL1 phosphorylation levels by LC-MS.\nThe qHtL1-GST and TaCDPK13-GST fusion proteins were obtained using a prokaryotic expression\nsystem, and then purified for in vitro phosphorylation assay. The qHtL1-GST fusion protein was mixed\nwith different amounts of TaCKPK13-GST, and phosphorylation levels of qHtL1 in each reaction were\nevaluated by mass spectrometry. In f1, no TaCDPK13-GST was added in the action system. From f2 to\nf4, the amounts of TaCDPK13-GST were gradually increased in the reaction system. By maintaining\nconsistent quantities of qHtL1-GST in each reaction, the amounts of four peptides containing identified\nphosphorylation sites were observed to increase in response to elevated levels of TaCKPK13-GST within\nthe reaction systems. Colored squares indicate varying degrees of phosphorylation, with dark blue\ndenoting low phosphorylation levels and deep red signifying high phosphorylation levels.\n4. Figure 5 illustrates the regulation of qHtL1 by TaOGT1; however, it is unclear how many\nreplicates were performed to ensure the reproducibility of the data. Furthermore, there\nseems to be an inconsistency in the interpretation of the data between Figure 5g and Figure\n5d. If Figure 5g shows a multi-fold increase in qHtL1-HA protein levels following\nvernalization treatment, how does this correlate with the consistent input levels shown in\nFigure 5d? This discrepancy raises questions about the experimental design and\nnormalization methods used. In addition, the number of replicates performed and the\nreproducibility of these results are not clearly stated, which needs to be addressed to\nincrease confidence in the results.\nResponse：\nThanks for your constructive comments. As you noted, different normalization methods\nwere used in the original Figure 5g and 5d. In the original Figure 5d, the amounts of qHtL1-\nHA under V0 and V21 conditions were normalized to be consistent to compare the levels\nof O-GlcNAc modified qHtL1-HA. In the original Figure 5g, the total protein levels\ncorresponding to V0, V7, V14 and V21 were maintained consistent and normalized to Actin,\nshowing a gradual increase of qHtL1-HA following extended vernalization durations.\nWe performed three biological replicates, which demonstrated that the levels of qHtL1-\nHA progressively increased with the duration of vernalization treatment. In addition, our\nresults from two biological replicates demonstrated that the O-GlcNAc modification\ndegrees of qHtL1-HA gradually increased following extended vernalization durations,\nwhich accompanied the elevation of total qHtL1-HA levels. The representative result is\npresented in Figure 7 below and in Figure 5d in the revised manuscript.\nConsidering consistent and partly overlapped results, the original Figure 5d and 5g\nwere removed in the revised manuscript.\nFigure 7. Vernalization induces accumulation of O-GlcNAcylated qHtL1 in vivo.\nChemoenzymatic labelling approach was conducted for analyzing O-GlcNAclated qHtL1 in wheat. The\nqHtL1-OE plants were vernalized at 4°C for 0, 7, 14, and 21 days, and then for total protein isolation and\nchemoenzymatic labelling. Constitutively expressed qHtL1 was detected using anti-HA antibody. The\nqHtL1 levels were normalized to Actin, and then for quantification of the percentage of O-GlcNAcylated\nqHtL1 relative its total amount. The relative protein abundance was quantified using ImageJ software.\n5. If Figure 5g shows that the qHtL1 level increases several fold between v0 and v21, it is\npuzzling why the activity of the aldolase level increases only very slightly (900 to 1200?).\nTo prove the model, is it possible to treat the wheat with the product (G3P) of qHtL1 and\nshow that it would affect the histone modification states on VRN1 and affect the flowering\ntime?\nResponse:\n(1) In enzymatic reactions, the augmentation of protein concentrations generally does not\ndemonstrate a linear correlation with the enhancement of enzyme activity. Despite a\nsignificant elevation in the protein level of qHtL1-HA relative to the V0 condition following\n21 days of vernalization, the enzyme activity did not exhibit a proportional increase. We\npostulate that this discrepancy is attributable to an insufficient rise in substrate levels.\nAccording to Figure 3f (original and revised versions), 21 days vernalization treatment\nresulted in a ~30% increase in aldolase activity in wheat qHtL1-OE plants compared to\nwild-type KN199, lead to early flowering phenotype. The degree of variation in qHtL1\nenzyme activity resulting in phenotypic modifications aligns with findings reported in a\nliterature, in which, a 33% reduction in aldolase activity significantly affects the intracellular\ngrowth of Salmonella (Ref # below).\n(2) We agree with your suggestion that, exogenous application of G3P, the product of\nqHtL1, should have an impact on histone modification states at VRN1 locus and flowering\ntime. It is reasonable to speculate that the substrate of qHtL1, FBP (fructose 1,6-\nbisphosphate), also possibly exert similar influence with G3P. We treated KN199 plants with\n25mM and 50mM FBP during vernalization period of V21, and found that FBP treated\nplants showed an early flowering time phenotype comparing with untreated control plants\n(Figure 8a, b below).\nIn addition, we further evaluate the role of FBP application on histone modification\nstates at VRN1 in KN199 wheat plants vernalized for 21 days (V21). The KN199 material,\nthree days post-germination, was subjected to a 21-day vernalization treatment. FBP was\napplied at 24 hours prior to sampling. H3K27ac and H3K14ac modification levels at region\nIV of the VRN-A1 locus were evaluated by ChIP-qPCR analysis, demonstrating that FBP\naddition resulted in a significant increase in these two regions (Figure 8c, d below).\nMoreover, quantitative PCR (qPCR) assays demonstrated a significant upregulation of the\ntranscriptional level of VRN1 in plants treated with FBP in comparison to the untreated\ncontrol (Figure 8e), indicating the role of H3K27ac and H3K14ac in activating VRN1.\nTogether these results support our proposed model in Figure 6f.\nRef #: Fitzsimmons, L., et al. (2018). \"Zinc-dependent substrate-level phosphorylation powers Salmonella\ngrowth under nitrosative stress of the innate host response.\" PLoS Pathog 14(10): e1007388.\n(a) Exogenous addition of FBP accelerates flowering in KN199. KN199 was treated with 25 mM and 50\nmM FBP during vernalization process, respectively. Wild-type KN199 plants and FBP-treated KN199\nplants were transplanted in a green house, following a 21-day vernalization treatment. Scale bar, 10 cm.\n(b) Heading time statistical analyses of (a). The box plots display the interquartile range, comprising the\nfirst quartile, median, and third quartile, while the whiskers extend from the minimum to the maximum\nvalues, (n=15 plants for each group). Two-tailed Student’s t-test for statistical analysis.\n(c) and (d) ChIP-qPCR assays of H3K14ac and H3K27ac levels in the region IV at VRN1. Schematic\ndiagram of VRN1, showing the region detected in ChIP-qPCR analysis. Two biological replicates were\nperformed, with each biological replicate comprising three technical replicates. Two-tailed Student’s t-test\nfor statistical analysis.\n(e) qRT-PCR assay shows that exogenous addition of FBP activates VRN1 transcription. The KN199\nplants, 3 days after germination, were subjected to a 21-day vernalization treatment. FBP was applied at\n24 hours prior to sampling. Data are means ± SD of three biological replicates. Two-tailed Student’s t-test\nwas used for statistical analysis.\n6. Additional Points:\nAll figures from Figure 1 through Figure 6 use traditional bar graphs. To improve visual\nclarity and allow for a better understanding of the actual data distribution, box plots or\nbeeswarm plots should be used.\nResponse：\nThank you for your valuable suggestion. In the revised manuscript, all bar graphs have\nbeen converted to box plots, with the exception of the qPCR results, UID mRNA-seq data\nand results of aldolase activity assay, which remain unchanged due to only three data\npoints.\nReviewer #1 (Remarks to the Author):\nHowever, I remain unconvinced that they have proven that\nTraesCS3A02G391100/FBA10/qHtL1 is the causal gene for the heading date phenotype.\nFor me, key issues are to do with the functional differences between the haplotypes and\nthe genetic confirmation that this is the causal gene. The genetic confirmation is\ncomplicated by the result they mention in the response document that a CRISPR-induced\nallele is lethal – so they cannot test its effect on heading date. I think that they should\ninclude this CRISPR allele and the lethality phenotype in the manuscript, so that the reader\nknows that this experiment was done. In the response document, they go on to say “In\ncontrast, two haplotypes of FBA10 (Hap1 vs Hap2) showed highly significant differences\nin heading date phenotypes (Figure 1c in the revised manuscript).” This difference seems\nstatistically supported. So, to be convincing they should show us that these two haplotypes\nhave polymorphisms in FBA10 consistent with the heading date phenotype. In other words,\nis the predicted amino acid sequence or the expression level of Hap2 of\nTraesCS3A02G391100/FBA10/qHtL1 consistent with altered activity of the gene and\ndelayed heading? I could not see these data anywhere, and without them I still find the\nmanuscript unconvincing.\nResponse:\nThank you for your comments and we greatly appreciate your hard work and constructive\nsuggestions, which greatly helped us to improve the manuscript.\n(1) We had previously attempted to simultaneously disrupt all A-, B- and D-homoeologs of\nFBA10/ qHtL1 in the KN199 background utilizing the CRISPR/Cas9 genome editing approach,\nand three single-guide RNAs (sgRNAs) were designed to target different sites in exon 2,\nrespectively, which was shown as Supplementary Figure 7 in the revised manuscript, as well\nas Response Figure1a and Response Figure2 below. PCR analysis utilizing primers specifically\ndesigned to amplify the Cas9 gene successfully identified four transgenic plants among the 32\nT0 generation plants (please see Response Figure1b below). Two independent mutant lines,\ncr-1 and cr-2 were identified by further sequencing analyses. Sequencing results revealed that\nthe cr-1 mutant exhibits a base substitution at 178 bp downstream of the ATG start codon within\nthe 3A genome, accompanied by a fragment deletion spanning from 181 bp to 248 bp. In the\n3B genome, a fragment deletion is detected between 180 bp and 246 bp downstream of the\nATG. In the 3D genome, a fragment deletion occurs between 177 bp and 330 bp downstream\nof the ATG. The cr-2 mutant demonstrated three distinct fragment deletions: one located\nbetween 181 bp and 248 bp downstream of the ATG start codon in the 3A genome, another\nspanning from 181 bp to 246 bp downstream of the ATG in the 3B genome, and a third ranging\nfrom 177 bp to 330 bp downstream of the ATG in the 3D genome. These data are presented in\nSupplementary Figure 7 of the revised manuscript, as well as Response Figure 1a and\nResponse Figure2 below.\nWe obtained only the two mutant lines previously mentioned, despite three transformation\nattempts had been conducted to generate FBA10/qHtL1 knockout mutants by CRISPR-Cas9\nsystem. The subsequent two attempts even failed to yield any positive transgenic plants. These\ntwo mutants exhibited markedly impaired growth and developmental phenotypes, ultimately\nresulting in lethal phenotype within two months. Unfortunately, we did not capture photographic\ndocumentation of the phenotypes of these mutants at that time. We suggest that the complete\nloss of FBA10/qHtL1 function results in a lethal phenotype. Consequently, we proceeded to\nconstruct RNAi materials targeting the FBA10/qHtL1 gene. The related information has been\nincorporated into the revised manuscript.\nResponse Figure 1 Generation and molecular identification of FBA10 knockout mutants\n(a) Three single guide RNA (sgRNA) sequences are designed to target regions within exon 2 of the three\nTraesCS3A02G391100/FBA10/qHtL1 homoeologs. In the diagram, exons are depicted as gray rectangles,\nand the positions of the target sites are indicated by red bold lines. Sequencing results highlight a base\nsubstitution in red and base deletions as red dashes, with numbers on the right indicating gap lengths.\n(b) PCR analysis to identify the presence of Cas9 in T0 plants produced via CRISPR-Cas9 genome editing.\nCas9-specific primers were designed to amplify a 338 bp fragment. A total of 32 T0 seedlings were\nobtained and screened, among which lines 1, 2, 16, and 30 were identified as positive transgenic plants.\nThese lines were subsequently subjected to sequencing analysis to confirm the mutation patterns.\nResponse Figure 2 Identification of cr-1 and cr-2 via sequencing\nNotes: In the diagram, exons are depicted as gray rectangles, and the positions of the target sequences\nare indicated by red bold lines. The target site sequences for the three sgRNAs and the mutation patterns\nin cr-1 and cr-2 mutants are shown. Sequencing results highlight a base substitution in red and base\ndeletions as red dashes, with numbers on the right indicating gap lengths.\n(2) We performed a comprehensive analysis of Hap-1 and Hap-2, revealing that most SNPs\nare situated within the promoter region (3A: 639082776-639084776) of the\nTraesCS3A02G391100/FBA10/qHtL1 gene. Utilizing the Binding Site Prediction tool available\non the PlantRegMap platform (https://plantregmap.gao-lab.org/index.php), we identified that\nthe SNP S3A_639082983 is located within the binding site (CACACACACACTCTATCTTTA)\nof the BBR-BPC transcription factor (Supplementary Data 4 in the revised manuscript). This\nbinding site has been previously validated to regulate gene expression in barley, as\ndocumented in the literature (Ref#1 below). Consequently, we propose that, despite the\nabsence of amino acid sequence alterations between Hap-1 and Hap-2, the SNP variation\nwithin the promoter region may influence the recognition of the BBR-BPC transcription factor,\nthereby modulating the expression level of FBA10/qHtL1. The information has been included\nin the revised manuscript.\nIt has been reported that FBA10/qHtL1 exhibits high expression levels at different\ndevelopment stages compared with its homologous family members, especially at germination\nand inflorescence emergence stages (Ref#2 below). We performed a new additional qRT-PCR\nanalysis to examine the differential expression of FBA10/qHtL1 in wheat plumules. This\nanalysis encompassed 12 accessions with the Hap-1 haplotype and 23 randomly selected\naccessions with the Hap-2 haplotype. The experiment was conducted with three independent\nbiological replicates. The qRT-PCR results indicated that the expression levels of FBA10/qHtL1\nare generally elevated in the Hap-1 group relative to the Hap-2 group, aligning with the\nphenotypic differences in heading time observed between the two haplotypes. The qRT-PCR\ndata were shown as Supplementary Figure 5 in the revised manuscript and Response Figure\n3 below. This result, in conjunction with other analyses, suggests that the predicted\nTraesCS3A02G391100/FBA10/qHtL1 is the candidate gene for the genomic region identified\nthrough GWAS analysis. The related information has been included in the revised version of\nthe manuscript.\nRef#1: Santi L., et al. (2003). The GA octodinucleotide repeat binding factor BBR participates in the\ntranscriptional regulation of the homeobox gene Bkn3. Plant Journal 34(6): 813-826.\nRef#2: Ref #1: Lv, G. Y., et al. (2017). Molecular characterization, gene evolution, and expression analysis\nof the fructose-1, 6-bisphosphate aldolase (FBA) gene family in Wheat (Triticum aestivum L.). Front Plant\nSci 8: 1030.\nResponse Figure 3 Expression analysis of FBA10/qHtL1 in Hap-1 and Hap-2 groups.\nNotes: The transcript levels of FBA10/qHtL1 in plumules of 12 Hap-1 and 23 randomly selected Hap-2\naccessions were analyzed by qRT–PCR. Three independent biological replicates were conducted, and\nthe data from all replicates are presented in the figure. Data are means ± s.d., two-tailed Student’s t-test\nwas used to determine the P value.\nResponse to Reviewer #1\nReviewer #1 (Remarks to the Author):\nIn this second revision, the authors have seriously considered my comments and edited the\nmanuscript accordingly. Including the data on the CRISPR mutants as well as the details of the\nsequence differences and expression levels of the haplotypes have improved the manuscript and\nanswered my criticisms. Showing that in Hap1 TraesCS3A02G391100/FBA10/qHtL1 is slightly\nmore highly expressed than in Hap2 in supplementary Figure 5 provides a logical basis for the\nheading date phenotype observed, although in my opinion the difference in mRNA level does look\nvery small. Nevertheless, I realize that the other two reviewers have already accepted the manuscript\nfor publication, and it does now present a logical argument supporting their conclusion that this is\nthe causal gene. Therefore, I support publication of this version.\nResponse:\nWe would like to express our sincere gratitude for your constructive feedback and for\nendorsing the acceptance of our manuscript. Your insights have been instrumental in\nhelping us improve our work to a deeper level."}, {"round": 2, "reviewer_comments": {"reviewer_1": "I appreciate that the authors have carefully considered my comments on the previous version and have provided a thoroughly argued response that includes new data. However, I remain unconvinced that they have proven that TraesCS3A02G391100/FBA10/qHtL1 is the causal gene for the heading date phenotype. For me, key issues are to do with the functional differences between the haplotypes and the genetic confirmation that this is the causal gene. The genetic confirmation is complicated by the result they mention in the response document that a CRISPR-induced allele is lethal – so they cannot test its effect on heading date. I think that they should include this CRISPR allele and the lethality phenotype in the manuscript, so that the reader knows that this experiment was done. In the response document, they go on to say “In contrast, two haplotypes of FBA10 (Hap1 vs Hap2) showed highly significant differences in heading date phenotypes (Figure 1c in the revised manuscript).” This difference seems statistically supported. So, to be convincing they should show us that these two haplotypes have polymorphisms in FBA10 consistent with the heading date phenotype. In other words, is the predicted amino acid sequence or the expression level of Hap2 of TraesCS3A02G391100/FBA10/qHtL1 consistent with altered activity of the gene and delayed heading? I could not see these data anywhere, and without them I still find the manuscript unconvincing.", "reviewer_2": "The authors have considered the comments and suggestions that I raised in the first round of review and have modified the manuscript accordingly.", "reviewer_3": "The authors responded satisfactory to my comments."}, "author_response": ""}, {"round": 3, "reviewer_comments": {"reviewer_1": "In this second revision, the authors have seriously considered my comments and edited the manuscript accordingly. Including the data on the CRISPR mutants as well as the details of the sequence differences and expression levels of the haplotypes have improved the manuscript and answered my criticisms. Showing that in Hap1 TraesCS3A02G391100/FBA10/qHtL1 is slightly more highly expressed than in Hap2 in supplementary Figure 5 provides a logical basis for the heading date phenotype observed, although in my opinion the difference in mRNA level does look very small. Nevertheless, I realize that the other two reviewers have already accepted the manuscript for publication, and it does now present a logical argument supporting their conclusion that this is the causal gene. Therefore, I support publication of this version."}, "author_response": ""}]

Decision letter