Synteny and distribution of genomic features of Q. glauca and L. edulis genomes.

(A) A phylogenetic tree of five Fagaceae species from three genera, Quercus, Lithocarpus, and Fagus. Bootstrap values for all branches are 100. (B) Synteny relationship between the Q. glauca and L. edulis genomes. The collinear blocks within the genomes were displayed by lines. Black lines show reverse alignment and gray lines represent regular. Numbers correspond to the 12 chromosomes. (C–D) Circos plots of the Q. glauca (C) and L. edulis genomes (D). Lanes depict circular representation of chromosome length (Mbp) (a), gene density (b), GC density (c), and repeat (d). Lines in the inner circle represent links between synteny-selected paralogs. (E) Gene copy-number variations between the Q. glauca and L. edulis genomes. The bar plot illustrates the number of genes in gene families categorized into 1:1 orthologs, N:N (many-to-many), 1:N (one-to-many), N:1 (many-to-one), and species-specific unique genes. (F–G) Mean gene expression levels in each copy number variation category for Q. glauca (F) and L. edulis (G). Different letters denote statistically significant differences (p < 0.05, Steel-Dwass test).

Seasonal gene expression dynamics.

(A) Hierarchical clustering of 213 gene expression profiles encompassing two tissues (leaf and bud) and four species (Q. glauca, Q. acuta, L. edulis, and L. glaber), collected over two-year. (B) Distribution of daily mean temperature on each sampling date across clusters, from L1 (top) to B3 (bottom). The Mann-Whitney U test was used to compare daily mean temperature between cluster L1 and L2 (p = 0.136). The Steel-Dwass test was applied for multiple comparisons across clusters B1–B3 (p < 0.05). Different letters denote statistically significant differences (p < 0.05). (C–D) Hierarchical clustering of 11749 genes based on the similarity of gene expression dynamics. Daily mean temperature, flowering and leaf flushing phenology (proportion) are shown in the upper three panels of the heatmaps.

Distribution of period and peak month in rhythmic gene expression.

(A– B) Distribution of period length (month) among the genes with significant rhythmicity calculated by RAIN for leaf (A) and bud (B) samples of each species. (C–D) Distribution of seasonal peaks among the genes with annual periodicity. The angle of each bar represents the peak month, while the bar height indicates the number of genes. Purple and green squares around each rose plot denote the observed phenology of flowering and leaf flushing (proportion), respectively.

Comparison of peak months in seasonal gene expression across species in buds.

(A) Intra-genus and (B) inter-genera comparisons of peak months in seasonal gene expression. The size and color of the circles represent the number of genes and the mean of Pearson’s correlation, respectively. Genes with negative correlations lower than −0.3 are highlighted in blue, while the month of leaf flushing is highlighted in green. (C) Plot of molecular phenology divergence index D across months. Triangles and circles indicate the values of D calculated from intra-genus and inter-genera comparison, respectively. The line and shaded envelope indicate the mean and standard deviation (s.d.). Different letters denote statistically significant differences (Nemenyi test, p < 0.05).

Phylogenetic constraints in the evolution of seasonal gene expression and the relationship between seasonal gene expression divergence and sequence evolution.

(A) Distribution of pairwise Pearson’s correlation coefficients for gene expression across species (Le: L. edulis, Lg: L. glaber, Qg: Q. glauca, Qa: Q. acuta; n = 11749). Colors indicate rhythmicity categories. Black bars show median correlation values, and different letters indicate significant differences (Nemenyi test, p < 0.05). (B) Heatmap of pairwise correlation coefficients for each gene, with mean gene counts in four seasonal expression categories: (1) annual rhythmic genes peaking in the growing season (March–November), (2) annual rhythmic genes peaking in winter (December–February), (3) half-annual rhythmic genes, and (4) arrhythmic genes. Genes are ordered by mean correlation across species pairs. (C–D) Conserved (C) and diverged (D) gene expression patterns, based on the top and bottom 5% of pairwise correlations. Solid lines and error bars indicate mean ± SD. Numbers in parentheses show peak month differences between Q. glauca and L. edulis. Peak month difference for BFT was not calculated due to arrhythmic expression in Q. glauca. (E) Relationship between evolutionary rate (dN/dS) and peak month difference (Δφ) for genes with annual periodicity in Q. glauca and L. edulis (n = 951). Black lines and transparent bands represent regression lines and 95% confidence intervals.