Multi-omics reveals divergent regulation of anthocyanin glycosylation and gibberellin biosynthesis underlying leaf variegation in Saxifraga stolonifera

Multi-omics reveals divergent regulation of anthocyanin glycosylation and gibberellin biosynthesis underlying leaf variegation in Saxifraga stolonifera | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Research Article Multi-omics reveals divergent regulation of anthocyanin glycosylation and gibberellin biosynthesis underlying leaf variegation in Saxifraga stolonifera Jianhang Zhang, jiecheng Li, Hai Xing, Feng Zhang, Chao Wang, and 3 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-9460163/v1 This work is licensed under a CC BY 4.0 License Status: Under Review Version 1 posted 8 You are reading this latest preprint version Abstract Saxifraga stolonifera Curtis is a shade-tolerant, variegated-leaf herb with ornamental and medicinal value that naturally displays three stable leaf phenotypes: green, white-variegated, and purple-variegated. These phenotypes simultaneously exhibit both air space-type and pigment-type variegations. However, the key pathways and regulatory networks underlying their formation remain largely unknown. Here, we generated a 2.01 Gb high-quality chromosome-scale genome assembly (2n = 36). Comparative genomic analyses revealed a recent whole-genome duplication specific to S. stolonifera following the shared γ-triplication event, along with 155 expanded gene families enriched in flavonoid and phenylpropanoid biosynthesis pathways. Metabolomic profiling identified 58 anthocyanin-related compounds, among which 10 key pigments including cyanidin-3-O-glucoside, peonidin-3-O-galactoside, and quercetin-3-O-glucoside, were responsible for the purple patches. Their accumulation corresponded with the up-regulation of the anthocyanidin 3-O-glucosyltransferase gene SsBZ1 ( Sst12G009310 , Sst14G006050 ). Hormone and transcriptome analyses showed that white-variegated leaves accumulate high levels of bioactive gibberellins GA4 and GA7, driven by increased expression of gibberellin 3β-dioxygenase SsGA3ox ( Sst05G017510 ). Weighted gene co-expression network analysis (WGCNA) further identified a GA-associated module (salmon) with SsGA3ox as the hub gene, promoting cell expansion and generating air spaces between epidermal and palisade tissues along the veins. Collectively, our high-quality genome, metabolome, and transcriptome resources demonstrate that SsGA3ox -mediated gibberellin biosynthesis drives air space-type leaf variegation, whereas SsBZ1 -controlled anthocyanin glycosylation produces pigment-type leaf variegation. These findings provide an integrative omics framework for dissecting leaf variegation mechanisms in Saxifraga and other ornamental plants. leaf variegation anatomy Saxifraga stolonifera gibberellins transcriptomics Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Figure 7 Introduction Variegated leaf plants represent a distinct group of colored-leaf species characterized by stable, heritable patterns formed by differentially pigmented leaf regions. These features make them valuable model systems for studying plant chromatology, mechanisms of adaptation to low-light environments, and chloroplast development[ 1 , 2 ]. These variegated leaves can be classified into five types: chlorophyll type, air space type, epidermis type, pigment type, and appendages type[ 3 ]. Among these, the air space type and pigment type are the most prevalent and often co-occur on the same leaf, in species such as Actinidia lanceolata , Erythronium dens-canis , and Saxifraga stolonifera [ 3 , 4 ]. Variegated leaves contribute to multiple adaptive functions, including acclimation to shade conditions, enhanced cold tolerance, defense against herbivores, and facilitation of self-pollination[ 2 , 5 – 8 ]. In addition, variegated plants possess significant ornamental value due to their visually striking foliage. Although naturally occurring leaf variegation has been studied for decades, previous research has largely focused on morphological classification and ecological functions; genomics, transcriptomics, and other omics studies have been limited by comparison. To date, transcriptomic studies of variegated leaves have primarily focused on the types of variegations associated with chlorophyll deficiency or red-leaf phenotypes. These investigations have elucidated the molecular mechanisms underlying chlorophyll degradation or anthocyanin accumulation in variegated tissues, as demonstrated in Ananas comosus var. bracteatus [ 9 ], Epipremnum aureum ‘Marble Queen’[ 10 ], Clivia miniata var. variegata [ 11 ], Alternanthera bettzickiana [ 12 ], Camellia sinensis ‘Zijuan’[ 13 ], and Perilla frutescens [ 14 ]. In transcriptome analyses of air space-type variegation in Trifolium pratense and Primulina pungentisepala , key pathways related to photosynthesis, redox regulation, cell-wall modification, and nitrogen metabolism have been examined[ 2 , 15 ]. Anthocyanins, a major class of water-soluble flavonoid pigments, play a critical role in determining the coloration of plant organs and tissues[ 16 , 18 ]. Under high-light conditions, the flavonoid biosynthesis pathway is activated, resulting in the accumulation of photoprotective and antioxidative flavonoids, particularly flavonols and anthocyanins, in leaf tissues[ 19 ]. The spatially restricted accumulation of anthocyanins is primarily responsible for pigment-type variegation[ 1 , 20 ]. For example, in Corydalis hemidicentra , insertion of a 254-bp transposon into the bHLH35 gene enhances anthocyanin biosynthesis, conferring an environmentally adaptive gray phenotype[ 21 ]. Despite these advances, the developmental dynamics and transcriptional regulatory networks governing both air space-type and pigment-type variegation remain poorly understood[ 2 ]. Saxifraga stolonifera is a shade-tolerant perennial exhibiting variegated leaves with both ornamental and medicinal properties, and it displays diverse leaf variegation patterns that support adaptation across variable environments[ 1 , 22 , 23 ]. In its natural habitat, this species commonly shows three stable leaf phenotypes: green-leaf type, white-variegated type, and purple-variegated type. Here, we report a high-quality genome assembly of S. stolonifera generated by integrating Illumina short-read sequencing, PacBio long-read sequencing, and high-throughput chromosome conformation capture (Hi-C) technologies. Furthermore, we performed comprehensive metabolomics and transcriptomic profiling of young leaves, first-stage adult leaves, and fully mature leaves across the three phenotypes. A weighted gene co-expression network analysis (WGCNA) was also conducted to identify key genes and regulatory modules associated with changes in hormone-mediated metabolites during leaf variegation. By integrating morphological and physiological traits with plant hormone profiles, metabolomics data, and transcriptomic landscapes, this study provides a high-quality reference genome and offers in-depth insights into the molecular mechanisms underlying leaf variegation. These resources establish a foundation for future genetic, genomic, and functional studies within the genus Saxifraga . Materials and methods Plant materials The green-leaf (G) and purple-variegated (P) phenotypes of S. stolonifera were collected from Tianmu Mountain, Zhejiang, China (30°18'30" N, 119°23'47″ E), while the white-variegated (W) phenotype was obtained from Doupeng Mountain, Guizhou, China (26°20′08″ N, 107°17′34″ E). The formal identification of S. stolonifera was performed by Dr. Jianhang Zhang. Voucher specimens were deposited in herbarium of East China Normal University (HSNU). As S. stolonifera is not a protected species, no collection permit was required. Starting in September 2020, all three phenotypes were cultivated in a greenhouse at the Biological Station of East China Normal University using a 1:3 mixture of nutrient soil and vermiculite as substrate. On April 20, 2021, healthy, uniformly growing individuals of each phenotype were selected for experimentation. From each plant, two leaves were sampled: one first-stage adult leaf and one fully developed adult leaf of similar size and growth status. The smallest leaf exhibiting typical morphological characteristics was designated as the first-stage adult leaf. On June, 22, 2025, the young leaves lacking characteristic pigmentation patterns were sampled. In this study, the second emerging leaf was defined as the young leaf, the third as the first-stage adult leaf (showing distinct patterning), and the sixth or seventh leaf as the mature adult leaf. Leaf developmental stages are illustrated in Fig. 1 . After harvest, leaves were wrapped in aluminum foil, rapidly frozen in liquid nitrogen, and stored at − 80°C for downstream analyses. In April 2023, individuals of all three phenotypes were transplanted from East China Normal University to the greenhouse at Shaoxing University. On May, 7, 2024, root tips were collected for chromosome karyotype analysis[ 24 ], and whole plants of the green-leaf phenotype at the flowering stage were sampled for genome sequencing. Macroscopic leaf morphology and phenotypic stability Eighteen plants, six representing each of the three leaf phenotypes with uniform growth characteristics, were selected for morphological analysis. To assess the adaptability of the green-leaf and purple-variegated types to varying light conditions, these individuals were exposed to distinct light intensities within the same greenhouse: full sunlight, shading with a single-layer shade net (reducing light intensity to approximately 25% of full sunlight), and shading with a double-layer shade net (reducing light intensity approximately 10%). Substrate composition, watering regimen frequency, and fertilization protocols were standardized across all experimental groups. Weekly visual documentation of leaf morphology and color variation was performed using a Canon EOS 800D camera, while leaf surface features were examined using a KEYENCE VHX-5000 stereoscope (Japan). Leaf structure and ultrastructure Fully developed leaves of the three S. stolonifera phenotypes exhibiting healthy and uniform growth were used for structural and ultrastructural analyses. Leaves were manually sectioned, and thin, intact slices were carefully selected to prepare temporary water-sealed slides. Cell morphology, dimensions, and arrangement; chloroplast number, distribution, and spatial organization; and the spatial arrangement of anthocyanin-containing cells were observed under bright-field mode at 4 ×, 10 ×, and 20 × magnification using an Echo RVL-100-M inverted integrated fluorescence microscope (Discover ECHO, US). Chloroplast autofluorescence was assessed in fluorescence mode via the CY5 channel to determine chloroplast abundance and intracellular localization. For scanning electron microscopy, leaf blades were longitudinally dissected along both sides of the midvein, preserving a central segment approximately 1-2cm wide to include both regular and variegated regions. Observations were conducted using a Hitachi S-4800 cold-field emission scanning electron microscope (Japan), operated at 1 kV. Measurements of the content of chlorophyll and anthocyanins The chlorophyll content was measured using a SPAD-502 Plus meter (Konica Minolta, Japan) and expressed as relative SPAD values. A total of 480 leaves were sampled from 40 uniformly growing plants representing the three leaf phenotypes. One leaf each from the first-stage adult and fully adult stages, matched for size and developmental status, was selected per plant. For fully adult leaves, measurements were taken either over the main vein (white-variegated area) or interveinal areas (purple-variegated area). Due to the difficulty of distinguishing variegated and non-variegated areas on first-stage adult leaves, measurements were restricted to the main vein area. Data were processed using Microsoft Excel (2019), and differences in the chlorophyll content among leaf phenotypes were analyzed by one-way ANOVA using R v4.1.2. The anthocyanin content was quantified using the MetWare platform (Wuhan, China; http://www.metware.cn/ ) and the AB Sciex QTRAP 6500 LC-MS/MS system. Nine S. stolonifera plants, three per phenotype, with similar growth were selected. One first-stage adult leaf and one fully adult leaf, matched for size and developmental stage, were collected from each plant, yielding nine samples in total. Samples were immediately frozen in liquid nitrogen and stored at -80 ℃ until analysis. The results of hierarchical clustering analysis (HCA) of both samples and metabolites were visualized as dendrograms accompanied by heatmaps. HCA was performed using the R package pheatmap with normalized metabolite signal intensities (unit variance scaling) represented on a color scale. Differentially abundant metabolites (DAMs) were identified based on pairwise comparisons between groups: adult leaves of the purple-variegated type (p) versus green-leaf type (g); white-variegated type (w) versus g; and p versus w. DAMs were considered significant if they met the following criteria: absolute Log 2 (FC) ≥ 1 (i.e., fold change ≥ 2 or≤ -2), P-value ≤ 0.05, and presence rate ≥ 1 across all samples. Identified metabolites were annotated using the KEGG COMPOUND database ( http://www.kegg.jp/kegg/compound/ ). Detection of the phytohormone content Phytohormone levels were quantified using an ultra-high-performance liquid chromatography-tandem mass spectrometry (UHPLC-MS/MS) system (ExionLC™ AD UHPLC-QTRAP 6500+, AB SCIEX Corp., Boston, MA, USA) at Novogene Co., Ltd. (Beijing, China). On 22 June 2025, young leaves and first-stage adult leaves were sampled. Nine healthy plants, three per phenotype (G, P, W), with consistent growth were selected. From each plant, one young leaf (second leaf) and one first-stage adult leaf (third leaf), comparable in size and developmental stage, were harvested, resulting in 18 samples. All samples were flash-frozen in liquid nitrogen and stored at -80℃ prior before analysis. HCA was performed on all detected plant hormones, with hormone concentrations normalized and clustered accordingly. Differentially abundant metabolites were identified through pairwise comparisons using the thresholds of fold change (FC) > 1.2 or FC < 0.833 and P-value < 0.05. Genome assembly and annotation Genome sequencing was conducted using the Illumina platform at Biomarker Technologies (Beijing, China), integrating PacBio HiFi, Hi-C, and RNA-seq data. For genome annotation, mixed RNA samples were collected from four healthy tissues (root, stem, leaf, and flower) and subjected to RNA sequencing. Genomic DNA was extracted from S. stolonifera leaves to construct 150-bp paired-end libraries, which were sequenced on the Illumina HiSeq X platform. Circular consensus sequencing (CCS) reads generated by PacBio HiFi were assembled using Hifiasm (version 0.19.9-r616; https://github.com/chhylp123/hifiasm ) with key parameters ‘-n 3 -I = 0’ and default settings. Chromosome-level scaffolding was performed using HapHiC ( https://github.com/zengxiaofei/HapHiC ) with default parameters. Genome completeness was evaluated using BUSCO ( https://github.com/metashot/busco ) in genome mode (-m genome). Jellyfish v2.3.05[ 25 ] was used to count k-mers (K = 31), and GCE v1.0.2[ 26 ] was applied to estimate genome size and heterozygosity from k-mer frequency distributions. Transposable elements (TEs) were identified using a combined homology-based and de novo approach. A de novo repeat library was constructed using RepeatModeler2 v2.0.1[ 27 ], which incorporates RECON v1.0.8[ 28 ] and RepeatScout v1.0.6[ 29 ]. Full-length LTR retrotransposons were initially detected using LTRharvest v1.5.10[ 30 ] and LTR_finder v1.07[ 31 ], and then refined with LTR_retriever. A non-redundant, species-specific TE library was generated by merging the de novo library with the Dfam v3.5 database. Final TE sequences were classified via RepeatMasker v4.12[ 32 ] using homology searches. Tandem repeats were annotated using TRF v4.09[ 33 ] and MISA v2.1[ 34 ]. Gene annotation was performed using the Makerpipeline[ 35 ]. To minimize annotation errors, gene predictions were performed using three complementary approaches: de novo prediction, transcript-based assembly, and homology-based alignment. Homology-based predictions utilized protein sequences from Arabidopsis thaliana , Bergenia scopulosa , Chrysosplenium sinicum , Liquidambar formosana , and Vitis vinifera , which were mapped to the S. stolonifera genome using GeMoMa v1.7[ 36 ]. De novo gene prediction was carried out using Augustus v3.2.3[ 37 ] and SNAP[ 38 ]. RNA-seq reads were aligned and assembled using HISAT2 v2.1.0[ 39 ] and StringTie v2.1.4[ 40 ], respectively, and used as evidence for gene prediction with GeneMarkS-T v5.1[ 41 ]. Additionally, PASA v2.4.1[ 42 ] was employed to predict gene models from RNA-Bloom v2.0.0[ 43 ] assembled unigenes and PacBio/ONT full-length transcripts. All predicted models were consolidated using EVM v1.1.1[ 44 ] and subsequently updated with PASA. The completeness of the annotated gene set was evaluated using BUSCO v5.2.2[ 45 ] in protein mode (-m proteins). Comparative genomic analysis A comparative genomic analysis was conducted involving S. stolonifera and 24 other plant species (Table S1 ). Protein sequences were clustered into gene families using OrthoFinder v2.4.0[ 46 ]. Functional annotations of gene families were derived from the Pfam V33.1 database[ 47 ]. Unique gene families in each species were identified through GO and KEGG enrichment analyses. Single-copy orthologous genes present in at least 80.0% of species (n = 799) were aligned using MAFFT v7.205 (--localpair --maxiterate 1000)[ 48 ]. A maximum-likelihood phylogenetic tree was inferred using iqtree v2.2.0 (JTT + F+I+G4, 1000 bootstrap replicates)[ 49 ], with Amborella trichopoda designated as the outgroup. Divergence times were estimated using the MCMCTree program in the PAML v4.9i package under default settings. Gene family expansions and contractions (family-wide P-values and Viterbi P-values < 0.05) were analyzed using CAFÉ v4.2.1[ 50 ]. The evolutionary dynamics of the S. stolonifera genome were further explored by calculating the synonymous substitution rate (Ks) for collinear gene pairs using WGDI v0.71[ 51 ]. RNA sequencing and bioinformatic analysis Total RNA was extracted from 38 samples, including one blind test sample from each of the P1 and W2 groups, comprising first-stage adult and fully adult leaves of the three leaf phenotypes. RNA purity and integrity were assessed using a NanoDrop ND-1000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA) and the Agilent RNA Nano 6000 Assay Kit on the Agilent Bioanalyzer 2100 (Agilent Technologies, Santa Clara, CA, USA). Qualified libraries were subjected to next-generation sequencing (NGS) on the Illumina platform at Biomarker Technologies (Beijing, China). Additional transcriptome data for young leaves, representing the three phenotypes before characteristic patterns become visible, were collected in 2025 and sequenced on the T7 platform at Novogene Co., Ltd. (Beijing, China). Clean reads were aligned to the S. stolonifera reference genome using HISAT2[ 39 ]. Read counts were quantified using the subread package in R (featureCounts), and gene expression levels were estimated using TPM (transcripts per million) via StringTie[ 40 ]. In the KEGG database, KO (KEGG Orthology) identifiers denote functionally conserved orthologous gene groups. Based on whole-genome functional annotation, all genes associated with the phenylpropanoid biosynthesis (ko00940), flavonoid biosynthesis (ko00941), and anthocyanin biosynthesis (ko00942) pathways were identified[ 16 , 52 – 55 ]. Genes involved in diterpenoid biosynthesis (map00904) were also annotated [ 56 ]. Protein domains of candidate genes were identified using HMMER via the Quick Gene Family Identification plugin in TBtools-II[ 57 ], and genes with inconsistent domain architectures were filtered. DEGs were identified from transcriptome data using DESeq2[ 58 ] with the following criteria: false discovery rate (FDR) < 0.05 and |log 2 fold change| ≥ 1. The Venny tool (accessed 27 October 2025; https://www.omicshare.com/tools/ ) was used to visualize overlapping DEGs. After discarding genes with low relative expression (FPKM > 1 in more than 90% of the samples), the WGCNA plugin (Format =normalized count, Normalized method = raw, Sample percentage = 0.9, Expression Cutoff = 1, Filter Method = MAD, Reserved genes Num. =20,000; R 2 coutoff = 0.8, Recommended = 6, Customized = 8; min Module Size = 30, module cuttree height = 0.25, select max blocksize = 20,000; x axis label angle = 45, KME cutoff = 0.2, Choose method = 2) in TBtools-II was used to identify modules of highly correlated genes attributable to plant hormones based on the filtered FPKM data[ 57 , 59 ]. Results Genome assembly and annotation of S. stolonifera To establish a genomic framework for dissecting leaf variegation mechanisms, we generated a high-quality chromosome-scale genome assembly of S. stolonifera. A green-phenotype individual of S. stolonifera from the Yonglai village population on Qingliangfeng Mountain in Jixi County, Anhui Province, China, was selected for sequencing (Fig. 2 a). Karyotype analysis confirmed diploidy with 2n = 36 chromosomes (Fig. 2 b). K-mer analysis estimated a genome size of ~ 1.63 Gb with 0.29% heterozygosity and 57.6% repetitive sequence, based on the 31-mer frequency distribution (Fig. S1 ). To achieve a high-quality genome assembly, we integrated multiple sequencing technologies by integrating Illumina short-read (116.53 Gb), PacBio HiFi long-read (99.65 Gb), and Hi-C (325.13 Gb) sequencing data (Table S1 ). The initial assembly yielded 2,139 contigs, which were subsequently anchored into 18 pseudo-chromosomes using Hi-C interaction maps, accounting for 94% of the assembly (Fig. 2 c, d; Tables 1 , S1, S2). The final genome assembly spans 2.01 Gb with a scaffold N50 of 3,535,068 bp and a BUSCO completeness score of 98.5%, indicating high contiguity and completeness (Table 1 ; Table S3). Table 1 Statistics for the final genome assembly of S. stolonifera Genome information S. stolonifera (PacBio + Hi-C) PacBio HiFi data (Gb) 99.66 Hi-C clean data (Gb) 325.13 Assembly size (Gb) 2.01 Number of contigs 2,139 Scaffold N50 (bp) 3,535,068 Scaffold max (bp) 14,862,665 Anchor ratio (%) 94.00% GC content 34.12% BUSCO (%) 98.50% Repetitive sequences constitute a major fraction of eukaryotic genomes and primarily include tandem repeats and interspersed repeats[ 60 , 61 ]. To annotate repetitive elements in the S. stolonifera genome, both de novo and homology-based prediction methods were employed. A total of 1,835,418 repetitive sequences were identified, spanning 822,048,413 bp (~ 73.49% of the assembled genome) (Table S4). Retroelements were the most abundant repeat class, accounting for 44.99% of the genome, followed by DNA transposons (14%) and tandem repeats (14.5%) (Table S4). Gene models for the S. stolonifera genome were predicted through an integrative approach combining ab initio, homology-based, and transcriptome-supported predictions. This strategy yielded a total of 37,191 protein-coding genes, with an average gene length of 3,914.91 bp and an average of 4.95 exons per gene (Fig. S2 ; Table S5). Of these, 36,664 genes were successfully annotated, resulting in an annotation rate of 98.58% (Tables S6, S7). TrEMBL analysis indicated that 36,585 genes (98.37%) were functionally annotated, with 28,614 encoding metabolic enzymes as classified by Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways (Table S7). The completeness and quality of the annotation were evaluated using BUSCO. Among the 1,614 conserved orthologs assessed, 1588 (98.39%) were complete, including 86.18% single-copy and 12.21% duplicated genes, indicating high annotation completeness (Table S8). This high-quality genome provides a robust reference for subsequent functional and evolutionary analyses. Genome evolutionary and whole-genome duplication analysis of S. stolonifera To place S. stolonifera in an evolutionary context, we compared its genome with 24 representative plant species (Table S9). A total of 4,130 shared gene families were identified across all species, while 402 gene families were unique to S. stolonifera , enriched in pathways related to aminoacyl-tRNA biosynthesis, riboflavin metabolism, phenylpropanoid biosynthesis, and phenylpropanoid / flavonoid biosynthesis (Fig. 3 a, b; Fig. S3; Tables S10, S11). These pathway enrichments suggest potential lineage-specific adaptations in secondary metabolite production in S. stolonifera . Phylogenomic analysis using 799 single-copy orthologs confirmed that Saxifragaceae is sister to Grossulariaceae, with S. stolonifera diverging from other Saxifragaceae species ( Bergenia scopulosa , Tiarella polyphylla , and Chrysosplenium sinicum ) approximately 59.39 million years ago (46.83–71.99 MYA) (Fig. 3 c). Gene family evolution analysis revealed 155 expanded and 17 contracted gene families in S. stolonifera (Fig. 3 c; Tables S12–S14). Notably, expanded families were significantly enriched in flavonoid and phenylpropanoid biosynthesis pathways (Fig. S4a; Table S13). Anthocyanin accumulation not only reduces photodamage[ 19 , 62 ] but also underlies purple pigmentation in plants[ 2 , 13 ]. Thus, the expansion of these gene families appears functionally linked to anthocyanin biosynthesis in S. stolonifera . Whole-genome duplications (WGDs) are recognized as key drivers of plant genome evolution[ 63 , 64 ]. Analysis of synonymous substitution rates (Ks) among paralogous gene pairs revealed two distinct peaks at ~ 0.32 and ~ 1.37, indicating two whole-genome duplication (WGD) events (Fig. 3 d; Fig. S5a), suggesting two independent WGD events. The ancient peak (Ks ≈ 1.37) corresponds to the core eudicot gamma whole-genome triplication, whereas the recent peak (Ks ≈ 0.32) represents a lineage-specific WGD that occurred after the divergence of S. stolonifera from other Saxifragaceae species (Fig. 3 e; Figs. S5b, S6)). Microsynteny analyses further support two rounds of WGD along the S. stolonifera lineage, with a 4:1 gene copy ratio observed between the S. stolonifera and Vitis vinifera genomes (Fig. 3 d, e). Higher chromosomal collinearity was observed between S. stolonifera and the three Saxifragaceae species ( T. polyphylla , B. scopulosa , and C. sinicum ), compared with R. nigrum (Grossulariaceae) (Fig. S5b), reinforcing the earlier divergence of S. stolonifera from R. nigrum and its closer relationship with the other Saxifragaceae members. Comparative genomic analyses between S. stolonifera and B. scopulosa or C. sinicum revealed syntenic depth ratios of 4:3 and 5:3, respectively (Fig. S6), indicating a shared WGD. Together, these findings provide strong evidence for the shared gamma-WGT event across Saxifragaceae and confirm an independent, lineage-specific WGD in S. stolonifera following this ancient duplication. Morphological and physiological basis of leaf variegation Transverse sections of the green-leaf type (G), white-variegated type (W), and purple-variegated type (P) leaves were prepared (Fig. 4 ; Fig. S7–S9). These sections revealed that air spaces were prevalent within the epidermal and palisade tissues along the central vein in both the white-variegated and purple-variegated types (Fig. 4 , b, c, e, f, h, i, k, l; Fig. S7–S8). Notably, in the white-variegated type, air spaces were observed among palisade cells (Fig. 4 b, e, h, k; Fig. S7). In contrast, epidermal and palisade tissue cells in the green-leaf type were densely packed (Fig. 4 a, d, g, j; Fig. S9). The chlorophyll content did not differ significantly among the three leaf phenotypes, nor across different regions of individual leaves or developmental stages (Fig. S10; Table S15). Therefore, the appearance of white veins primarily results from air spaces located between the epidermal and palisade tissues along the central vein in the variegated types. In transverse sections, the spatial distribution of red pigments varied across leaf phenotypes. In the white-variegated type, pigment-containing cells were scattered throughout the mesophyll and hypoepidermal tissues (Fig. S7c), whereas in the green-leaf type, they were predominantly confined to the hypoepidermis (Fig. S9c). Neither the green-leaf nor the white-variegated type exhibited purple patches on the adaxial surface (Fig. 4 a–c; Fig. S7a, S9a). Conversely, the purple-variegated type displayed accumulation of red pigments in the lower layers of the palisade tissue, leading to distinct purple patches visible on the upper leaf surface (Fig. S8a–e). Collectively, the prominent white veins in the white-variegated type arise from intercellular air spaces between the epidermis and palisade tissue, coupled with altered cellular organization of the palisade layer. The conspicuous purple patches in the purple-variegated type result from the localized accumulation of red pigments in sub-palisade cells, while the associated white venation shares the same structural origin as in the white-variegated type. According to the classification of leaf variegations by Zhang et al .[ 1 ], the white veins in the white-variegated type represent an air space-type variegation, whereas the purple-variegated type combines both pigment-based and air space-type mechanisms. Metabolomic profiling identifies anthocyanin glycosides as key pigments for purple-variegated type leaves To identify the metabolites responsible for purple-variegated type, we profiled anthocyanins in leaves of the three phenotypes using LC-MS/MS. A total of 58 anthocyanin-related compounds were detected (Table S16). Comparative analysis revealed 10 differentially accumulated metabolites (DAMs) enriched in the anthocyanin and flavonol biosynthesis pathways (ko00942, ko00944), including cyanidin-3-O-glucoside, peonidin-3-O-galactoside, and quercetin-3-O-glucoside (Fig. 5 a; Table S17). Notably, these compounds accumulated at significantly higher levels in purple-variegated type compared to green-leaf type, white-variegated type leaves (Fig. 5 a; Table S17). All identified pigments are glycosylated forms of anthocyanidins, suggesting that glycosylation plays a critical role in their stabilization and accumulation. These results pinpoint anthocyanin glycosides as the chemical basis of the purple patches. Plant hormone profiling identifies for white-variegated type leaves Plant hormones profiling detected 46 hormones across eight classes (Fig. 5 b; Table S18). Comparative analysis identified 35 differentially accumulated hormones, with the most pronounced differences observed between purple-variegated type and green-leaf type leaves (Table S19). Strikingly, both purple-variegated type and white-variegated type (W and P) accumulated significantly higher levels of multiple gibberellins (GA4, GA7, GA9, GA15, GA24) compared to green-leaf type leaves, in both young and first-stage adult leaves (Fig. 5 c; Table S19). Among these, GA4 and GA7 are bioactive forms known to promote cell expansion[ 65 ]. The elevated GA levels in variegated leaves point to a potential role for gibberellin in mediating air space formation. Transcriptomic analysis identifies SsBZ1 and SsGA3ox as key regulators of pigment type and air space type leaf variegation To identify genes associated with the formation of air space-type and pigment-type variegation, RNA sequencing was performed on leaves across the three developmental stages of the three leaf phenotypes in S. stolonifera . An average of 7.32 Gb clean reads per sample was obtained (Table S20), with an average mapping rate of 85%. A total of 18,923 differentially expressed genes (DEGs) were identified across all comparisons, with the highest number observed in the W2 vs. G2 comparison (6,209 DEGs; Fig. S11a; Table S21). Among these, 4,457 DEGs were commonly differentially expressed in at least two of the three stage specific comparisons between the purple-variegated and green-leaf types (YP-vs-YG, P1-vs-G1, and P2-vs-G2; Fig. S11b; Table S21). Similarly, 4,895 DEGs were shared across at least two comparisons between the white-variegated and green-leaf types (YW-vs-YG, W1-vs-G1, and W2-vs-G2; Fig. S11d; Table S21). Integrating these results with our high-quality S. stolonifera genome, we focused on key regulatory genes involved in anthocyanin and gibberellin biosynthesis, reconstructing the anthocyanin (Fig. 6 ) and gibberellin (Fig. 7 a) pathways associated with leaf variegation. We annotated 63 enzymatic genes involved in the anthocyanin biosynthesis pathway and visualized their expression patterns across three developmental stages of the three leaf phenotypes (Table S22). By reconstructing anthocyanin synthesis in S. stolonifera leaves, we found that Ss4CL ( Sst04G018150 , Sst04G018170 , Sst04G019880 , Sst14G006060 ), SsCHS ( Sst12G013330, Sst08G008480 ), and SsBZ1 ( Sst14G006050, Sst12G009310 ) were expressed at relatively high levels in fully developed leaves of the purple-variegated type (P2) (Fig. 6 ; Fig. S11 c; Table S22, S23). Notably, anthocyanidin 3-O-glucosyltransferase genes ( SsBZ1 : Sst12G009310 , Sst14G006050 ) showed elevated expression in both the white-variegated (W2) and purple-variegated (P2) leaf types, coinciding with the significant accumulation of key pigmentation metabolites: cyanidin-3,5-O-diglucoside, cyanidin-3-O-glucoside, peonidin-3-O-galactoside, peonidin-3-O-glucoside, peonidin-3,5-O-diglucoside, and quercetin-3-O-glucoside (Fig. 5 a; Table S17). In gibberellin biosynthesis, 24 enzymatically active genes were identified, and their expression profiles were analyzed across the three leaf phenotypes (Table S24). Reconstruction of the gibberellin pathway in S. stolonifera revealed that gibberellin 3β-dioxygenase ( SsGA3ox : Sst05G017510 ) was highly expressed in all three phenotypes, particularly in the white-variegated type (groups W1 and W2), likely contributing to sustained high concentrations of bioactive GA4 (Fig. 7 a; Fig. S11 e; Table S25). These expression patterns are consistent with the significant accumulation of GA4, the primary active gibberellin, in these tissues (Fig. 5 b, c, 7 a; Table S19). Collectively, these results suggest that the upregulated expression of GA 3-oxidases ( GA3ox ), which maintains elevated levels of active gibberellin, is a key factor underlying white vein formation. WGCNA identifies a gibberellin-associated co-expression module with SsGA3ox To gain further insights into the regulatory mechanisms underlying gibberellin content fluctuations during leaf variegation, weighted gene co-expression network analysis (WGCNA) was performed. A total of 9,960 DEGs were filtered and retained based on Fragments Per Kilobase of transcript per Million mapped reads (FPKM) > 1 in more than 90% of the samples and showing significant differential expression across the three developmental stages of the three leaf phenotypes (Fig. S12). This analysis revealed 23 co-expression modules (Fig. 7 b, c, color-coded). Among these gene co-expression subnetworks (Fig. 7 c), several exhibited strong correlations with the synthesis and activation of the five gibberellins. Notably, the two bioactive gibberellins, GA4 and GA7, showed highly significant positive correlations with the salmon-coded module, with correlation coefficients of 0.9 and 0.94, respectively. Gene expression patterns within this module were characterized by high transcript levels in the white-variegated and purple-variegated leaf types (Fig. 7 a; Fig. S11 d, e). Within the salmon-coded module, the key gene SsGA3ox ( Sst05G017510 ) displayed significantly elevated expression in white-variegated leaves, particularly in fully expanded adult leaves (Fig. 7 a; Fig. S13). These findings support the hypothesis that upregulated SsGA3ox expression maintains elevated levels of active gibberellins, thereby promoting rapid leaf expansion and contributing to the formation of intercellular air spaces between the epidermal and palisade tissues along the veins in air space type-leaf variegation. Discussion Here, we present a high-quality chromosome-scale genome assembly for S. stolonifera . Ks distribution analyses provided compelling evidence for a shared gamma WGT event among Saxifragaceae species and revealed an independent WGD event in S. stolonifera after the gamma-WGT. Comparative analysis of gene families across 25 species identified expanded gene families in S. stolonifera that are highly enriched in phenylpropanoid biosynthesis, fatty acid metabolism, stilbenoid diarylheptanoid and gingerol biosynthesis, and flavonoid biosynthesis. Additionally, we found that GA3ox genes in S. stolonifera play a significant role in maintaining elevated levels of active gibberellins. Collectively, these findings suggest that high concentrations of active gibberellins may have been crucial for the development of air spaces between the epidermal and palisade tissues along leaf veins. Localized anthocyanin accumulation determines the pigment type of S. stolonifera Bioinformatic analysis of S. stolonifera leaves indicated that the up-regulated expression of the anthocyanidin 3-O-glucosyltransferase gene ( SsBZ1 : Sst12G009310 , Sst14G006050 ) was a primary contributor to the formation of the purple-variegated phenotype. Anthocyanidin 3-O-glucosyltransferase is a key glycosyltransferase in plants[ 66 ]. Glycosylation plays a critical role in secondary metabolite biosynthesis by enhancing the stability, solubility, subcellular localization, and biological activity of conjugated compounds[ 66 – 69 ]. Enhanced expression of this gene promotes the glucosylation of anthocyanidins into stable anthocyanins, thereby preserving their coloration[ 55 ]. Previous studies have shown that the anthocyanidin 3-O-glucosyltransferase gene is up-regulated in anthocyanin-rich, non-green leaves, such as those of Alternanthera bettzickiana [ 12 ], Camellia sinensis ‘Zijuan’[ 13 ], and Perilla frutescens [ 14 ]. Similarly, in transcriptomic analyses of pigment-type variegation in Begonia masoniana , the UFGT gene (anthocyanidin 3-O-glucosyltransferase), involved in anthocyanin glycosylation is also up-regulated[ 17 ]. In this study, expression of SsBZ1 ( Sst12G009310 , Sst14G006050 ) was significantly up-regulated in the purple-variegated type compared with the green-leaf type (Fig. 6 ; Table S22). The content of visible anthocyanins was markedly higher in the purple-variegated type than in the green-leaf type, consistent with the elevated expression of SsBZ1 (Fig. 5 a, 6 ; Fig. S11 c; Table S17). Therefore, up-regulation of anthocyanidin 3-O-glucosyltransferase likely underlies the formation of the purple-variegated phenotype. Furthermore, compared with the white-variegated type, visible anthocyanins in the purple-variegated type were localized specifically within palisade tissue cells between the veins rather than scattered in the mesophyll cells (Fig. 4 ; Fig. S7–9). This spatial pattern explains why both types contain similar levels of visible anthocyanins (Fig. 5 a), yet only the purple-variegated type exhibits distinct purple patches. The localized accumulation of anthocyanins in palisade cells constitutes the structural basis for the observed pigmentation pattern and has also been documented in other pigment-variegated plants, including Nephelaphyllum pulchrum , Streptolirion volubile , and Vriesea saundersii [ 1 ]. In conclusion, up-regulation of the anthocyanidin 3-O-glucosyltransferase gene enhances anthocyanin accumulation, forming the molecular foundation of pigment-type variegation in S. stolonifera . However, the molecular mechanisms governing the specific localization of anthocyanins in palisade tissue cells between veins remain unclear. High gibberellin levels promote air spaces formation Based on the chromosome-level genome sequence, combined with gene co-expression and plant hormone profiling, we reconstructed the gibberellin biosynthesis pathway and its regulatory network in S. stolonifera . GAs are phytohormones that regulate multiple aspects of plant development[ 70 ], including leaf expansion[ 71 ], morphogenesis[ 72 ], and final leaf size[ 73 ]. Elevated GA concentrations accelerate cell division and leaf elongation rates[ 74 – 76 ]. Hormone profiling revealed that increased levels of bioactive GAs, particularly GA4 and GA7, play a pivotal role in the formation of air space-type leaf variegation in S. stolonifera (Fig. 5 c). During early leaf development, significant accumulation of multiple gibberellins, including GA4, GA7, GA9, GA15, and GA24, was observed in both the purple-variegated and white-variegated types (Fig. 5 c; Table S19). GA3ox catalyzes the final step in gibberellin biosynthesis, converting precursor GAs into the bioactive forms GA1 and GA4[ 77 – 82 ]. These enzymes are also known as gibberellin 3β-dioxygenases. Transcriptomic data demonstrated that the gibberellin 3β-dioxygenase gene ( SsGA3ox: Sst05G017510 ), which mediates this conversion, was highly expressed in white-variegated leaves, especially at the fully mature stage (Fig. 7 a; Fig. S11d e, S13). This sustained expression likely maintains elevated levels of active gibberellins, promoting rapid cell division and leaf expansion, key processes facilitating physical separation between epidermal and palisade tissue cells and resulting in air space formation along veins. Supporting this, WGCNA identified a co-expression module (salmon-colored) strongly correlated with GA4 and GA7 levels (Fig. 7 c), and genes within this module, including SsGA3ox , were highly expressed in variegated leaves (Fig. S12, S13). Together, these results indicate that gibberellin mediated cell expansion and tissue patterning form the structural basis of air space-type variegation. This mechanism enhances our understanding of leaf morphological diversity and provides a physiological and molecular explanation for the development of non-pigmented white veins in variegated leaves of S. stolonifera . In sum, the high-quality reference genome of S. stolonifera presented here serves as a valuable resource for studying evolutionary dynamics and phenotypic diversity within the genus Saxifraga . Moreover, comprehensive transcriptomic, anthocyanin metabolic, and plant hormone profiling analyses provide novel insights into the molecular mechanisms underlying both pigment-type and air space-type variegation. The data generated in this study represent essential resources for future functional genomics and genetic investigations in Saxifraga . Declarations Ethics approval and consent to participate Not applicable Consent for publication All authors agree to the publication of this manuscript upon its acceptance. Conflict of interest The authors declare no competing interests. Availability of data and materials The data that support the findings of this study have been deposited in the CNSA (accessed on 28 November 2024, https://db.cngb.org/cnsa/) of CNGBdb with accession code CNP0006535. Acknowledgments We extend our sincere thanks to Dr. Biao Xiong (Guizhou University) for guiding some data analyses. Our gratitude goes to Quangang Xue, Jiayu Jin, and Boyang Xie (Shaoxing University) for their kind help with the materials planting. We thank TopEdit (www.topeditsci.com) for linguistic assistance during the preparation of this manuscript. We also appreciate the assistance of Dr. Bin Dong (Zhejiang A&F University) for his help in improving our manuscript. AI declaration During the preparation of this work, the author(s) used [DeepSeek] solely for language polishing and to improve readability. After using this tool, the author(s) reviewed and edited the content as needed and take(s) full responsibility for the content of the publication. Author contributions H.L., S.L., J.L., and J.Z. conceived and designed the project; H.L., J.Z., and S.L. provided financial support; S.L., J.L., and J.Z. collected the samples; S.L., J.L., and J.Z. performed the experiments; F.Z., H.X., F.T., J.Z., and S.L. analyzed the data; F.Z., H.X., F.T., and C.W. contributed to project discussion. H.L., J.Z., and S.L. wrote the manuscript. All authors contributed critically to the drafts and gave final approval for publication. All authors have read and agreed to the published version of the manuscript. Funding The project was supported by a grant from the National Natural Science Foundation of China (32300178) to S.L., the Startup Fund for Shaoxing University (13011001002/218) to J.Z., Biodiversity Investigation, Observation and Assessment Program (2019-2023) of Ministry of Ecology and Environment of China to H.L. Supplemental data Supplementary data are available at BMC Genomics online. References Zhang JH, Zeng JC, Wang XM, et al. A revised classification of leaf variegation types. Flora. 2020;272:151703. Zhang JH, Li JC, Zou L, et al. Transcriptome analysis of air space-type variegation formation in Trifolium pratense . Int J Mol Sci. 2022;23:7794–7794. Zhang JH. Classification of leaf variegation types in angiosperm and variegation related genes of Trifolium pretense . East China Normal Univ. 2022. La Rocca N, Pupillo P, Puppi G et al. Erythronium dens-canis L. (Liliaceae): An unusual case of change of leaf mottling. Plant Physiol Biochem. 2014; 74: 108–117. Smith AP. Ecology of a leaf color polymorphism in a tropical forest species: habitat segregation and herbivory. Oecologia. 1986;69:283–7. Sheue CR, Pao SH, Chien LF, et al. Natural foliar variegation without costs? The case of Begonia . Ann Bot. 2012;109:1065–74. Pao SH, Chesson P, Peng CI et al. Characteristics and ecological significance of natural foliar variegation from Begonia formosana . In: The 99th Annual Meeting, 10th-15th (2014) August. Sacramento, California, America. Ecological Society of America. Shelef O, Summerfield L, Lev-Yadun S, et al. Thermal benefits from white variegation of Silybum marianum leaves. Front Plant Sci. 2019;10:688. Li X, Kanakala S, He YH, et al. Physiological characterization and comparative transcriptome analysis of white and green leaves of Ananas comosus var. bracteatus . PLoS ONE. 2017;12:e0169838. Sun YH, Hung CY, Qiu J, et al. Accumulation of high OPDA level correlates with reduced ROS and elevated GSH benefiting white cell survival in variegated leaves. Sci Rep. 2017;7:1–16. Wang QM, Cui JG, Dai HY et al. Comparative transcriptome profiling of genes and pathways involved in leaf-patterning of Clivia miniata var. variegata . Gene . 2018; 677: 280–8. Li WJ, Li HG, Shi LS, et al. Leaf color formation mechanisms in Alternanthera bettzickiana elucidated by metabolite and transcriptome analyses. Planta. 2022;255:59. Xie H, Zhu JY, Wang H, et al. An enhancer-transposable element from purple leaf tea varieties underlies the transition from evergreen to purple leaf color. Plant Commun. 2025;6(2):101176. Liu XN, Zhai YN, Liu JY, et al. Comparative transcriptome sequencing analysis to postulate the scheme of regulated leaf coloration in Perilla frutescens . Plant Mol Biol. 2023;112:119–42. Chen J, Li Y, He D, et al. Cytological, physiological and transcriptomic analysis of variegated Leaves in Primulina pungentisepala offspring. BMC Plant Bio. 2022;22:419. Tanaka Y, Sasaki N, Ohmiya A. Biosynthesis of plant pigments: anthocyanins, betalains and carotenoids. Plant J. 2008;54:733–49. Yang T, Xue ZZ, Li N. Reference genes selection and validation in Begonia masoniana leaves of different developmental stages. Acta Hortic Sin. 2021;48:2251–61. Zhu J, Wang YZ, Wang QY, et al. The combination of DNA methylation and positive regulation of anthocyanin biosynthesis by MYB and bHLH transcription factors contributes to the petal blotch formation in Xibei tree peony. Hortic Res. 2023;10:uhad100. Araguirang GE, Richter AS. Activation of anthocyanin biosynthesis in high light-what is the initial signal? New Phytol. 2022;236:2037–43. Cui WH, Guan KY. Diversity of leaf variegation in Chinese Begonias. Plant Divers Resour. 2013;35:119–27. Zhang H, Zhang P, Niu Y, et al. Genetic basis of camouflage in an alpine plant and its long-term co-evolution with an insect herbivore. Nat Ecol Evol. 2025;9:628–38. Kong LY, Ou GJ, Wang CH. Application of Saxifraga stolonifera Curt. to landscape architecture. Shaanxi For Sci Technol. 2014;2:77–8. Tang SM, Cai WQ, Zhang DM, et al. Effects of light intensites on morphological and physiological idexes of three ornamental cultivars of Saxifraga stolonifera . Guihaia. 2023;43:699–711. Li GF, He Y, Qin SH, et al. Optimization of chromosome preparation and karyotype analysis of polyploid watermelons. J South Agric. 2017;48:663–8. Marçais G, Kingsford C. A fast, lock-free approach for efficient parallel counting of occurrences of k-mers. Bioinformatics. 2011;27:764–70. Liu B, Shi Y, Yuan J, et al. Estimation of genomic characteristics by analyzing k-mer frequency in de novo genome projects. arXiv. 2012. 10.48550/arXiv.1308.2012.1308.2012v2 . Flynn JM, Hubley R, Goubert C, et al. RepeatModeler2 for automated genomic discovery of transposable element families. PNAS. 2020;117:9451–7. Bao Z, Eddy SR. Automated de novo identification of repeat sequence families in sequenced genomes. Genome Res. 2002;12:1269–76. Price AL, Jones NC, Pevzner PA. De novo identification of repeat families in large genomes. Bioinformatics. 2005; Suppl 1: i351. Ellinghaus D, Kurtz S, Willhoeft U. LTRharvest, an efficient and flexible software for de novo detection of LTR retrotransposons. BMC Bioinforma. 2008;9:18. Xu Z, Wang H. LTR_FINDER: an efficient tool for the prediction of full-length LTR retrotransposons. Nucleic Acids Res. 2007;35:W265–8. Tarailo-Graovac M, Chen NS. Using RepeatMasker to identify repetitive elements in genomic sequences. Curr Protocols Bioinf 2009; 4: 4.10. 11–14.10. 14. Benson G. Tandem repeats finder: a program to analyze DNA sequences. Nucleic Acids Res. 1999;27:573–80. Beier S, Thiel T, Münch T, et al. Misa-web: a web server for microsatellite prediction. Bioinformatics. 2017;33:2583–5. Cantarel BL, Korf I, Robb SMC. MAKER: An easy-to-use annotation pipeline designed for emerging model organism genomes. Genome Res. 2008;18:188–96. Keilwagen J, Wenk M, Erickson JL, et al. Using intron position conservation for homology-based gene prediction. Nucleic acids res. 2016;44:e89–89. Stanke M, Keller O, Gunduz I, et al. AUGUSTUS: ab initio prediction of alternative transcripts. Nucleic Acids Res. 2006;34:W435–9. Korf I. Gene finding in novel genomes. BMC Bioinforma. 2004;5:59. Kim D, Langmead B, Salzberg SL. HISAT: A fast spliced aligner with low memory requirements. Nat Methods. 2015;12:357–60. Pertea M, Kim D, Pertea GM, et al. Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat Protoc. 2016;11:1650–67. Tang S, Lomsadze A, Borodovsky M. Identification of protein coding regions in RNA transcripts. Nucleic Acids Res. 2015;43:e78. Haas BJ, Delcher AL, Mount SM, et al. Improving the Arabidopsis genome annotation using maximal transcript alignment assemblies. Nucleic Acids Res. 2003;31:5654–66. Nip KM, Chiu R, Yang C, et al. RNA-Bloom enables reference-free and reference-guided sequence assembly for single-cell transcriptomes. Genome Res. 2020;30:1191–200. Haas BJ, Salzberg SL, Zhu W, et al. Automated eukaryotic gene structure annotation using EVidenceModeler and the Program to Assemble Spliced Alignments. Genome Biol. 2008;9:R7. Simão FA, Waterhouse RM, Ioannidis P, et al. BUSCO: assessing genome assembly and annotation completeness with single-copy orthologs. Bioinformatics. 2015;31:3210–2. Emms DM, Kelly S. OrthoFinder: phylogenetic orthology inference for comparative genomics. Genome boil. 2019;20:1–14. Mistry J, Chuguransky S, Williams L, et al. Pfam: The protein families database in 2021. Nucleic Acids Res. 2021;49:D412–9. Katoh K, Asimenos G, Toh H. Multiple alignment of DNA sequences with MAFFT. Methods Mol Biol. 2009;537:39–64. Nguyen LT, Schmidt HA, von Haeseler A, et al. IQ-TREE: a fast and effective stochastic algorithm for estimating maximum-likelihood phylogenies. Mol Biol Evol. 2015;32:268–74. Han MV, Thomas GWC, Lugo-Martinez J, et al. Estimating gene gain and loss rates in the presence of error in genome assembly and annotation using CAFE 3. Mol Biol Evol. 2013;30:1987–97. Sun PC, Jiao BB, Yang YZ, et al. A user-friendly toolkit for evolutionary analyses of whole-genome duplications and ancestral karyotypes. Mol Plant. 2022;15:1841–51. Meng Y, Zhang H, Fan Y, et al. Anthocyanins accumulation analysis of correlated genes by metabolome and transcriptome in green and purple peppers ( Capsicum annuum ). BMC Plant Biol. 2022;22:358. Wang Q, Zhu J, Li B, et al. Functional identification of anthocyanin glucosyltransferase genes: a Ps3GT catalyzes pelargonidin to pelargonidin 3-O-glucoside painting the vivid red flower color of Paeonia. Planta. 2023;257:56. Geng Z, Liu M, Wang Y, et al. Transcriptomic and metabolomic analyses reveal CmMYB308 as a key regulator in the pink flower color variation of ‘Dante Purple’ chrysanthemum. Plant Cell Rep. 2024;43:157. Xia XM, Du HL, Hu XD, et al. Genomic insights into adaptive evolution of the species-rich cosmopolitan plant genus Rhododendron . Cell Rep. 2024;43:114745. Liu FH, Zhao Y, Wang XR, et al. Regulatory roles of GA and sugar metabolism genes with metabolomic profiling during seed dormancy release in Gleditsia sinensis . BMC Plant Biol. 2025;25:1231. Chen C, Wu Y, Li J, et al. TBtools-II: a one for all, all for one bioinformatics platform for biological big-data mining. Mol Plant. 2023;16:1733–42. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550. Wang X, Liang S, Yang WQ, et al. MetMiner: A user-friendly pipeline for large-scale plant metabolomics data analysis. J Integr Plant Biol. 2024;66:2329–45. Bendich AJ, Rogers SO. The biological and evolutionary consequences of competition between DNA sequences that benefit the cell and DNA sequences that benefit themselves. Nucleic Acids Res. 2025;53:gkaf589. Wang F, Jiang Z, Gao J, et al. Analysis of the Rehmannia chingii geneome identifies RcCYP72H7 as an epoxidase in iridoid glycoside biosynthesis. Nat Commun. 2025;16:6035. Solovchenko AE, Merzlyak MN. Screening of visible and UV radiation as a photoprotective mechanism in plants. Russ J Plant Physiol. 2008;55:719–37. Jiao YN, Wickett NJ, Ayyampalayam S, et al. Ancestral polyploidy in seed plants and angiosperms. Nature. 2011;473:97–100. Adams K. Genomic clues to the ancestral flowering plant. Science. 2013;342:1456–7. Hauvermale AL, Ariizumi T, Steber CM. Gibberellin signaling: a theme and variations on DELLA repression. Plant physiol. 2012;160:83–92. Yin RH, Messner B, Faus-Kessler T, et al. Feedback inhibition of the general phenylpropanoid and flavonol biosynthetic pathways upon a compromised flavonol-3-O-glycosylation. J Exp Bot. 2012;63:2465–78. Li Y, Baldauf S, Lim EK, et al. Phylogenetic analysis of the UDP-glycosyltransferase multigene family of Arabidopsis thaliana . J Biol Chem. 2001;276:4338–43. Meßner B, Thulke O, Schäffner AR. Arabidopsis glucosyltransferases with activities toward both endogenous and xenobiotic substrates. Planta. 2003;217:138–46. Bowles D, Lim EK, Poppenberger B, et al. Glycosyltransferases of lipophilic small molecules. Annu Rev Plant Biol. 2006;57:567–97. Salazar-Cerezo S, Martínez-Montiel N, García-Sánchez J, et al. Gibberellin biosynthesis and metabolism: A convergent route for plants, fungi and bacteria. Microbiol Res. 2018;208:85–98. Li J, Gao HW, Jiang JS, et al. Overexpression of the Galega orientalis gibberellin receptor improves biomass production in transgenic tobacco. Plant Physiol Biochem. 2013;73:1–6. Gao LW, Lyu SW, Tang J, et al. Genome-wide analysis of auxin transport genes identifies the hormone responsive patterns associated with leafy head formation in Chinese cabbage. Sci Rep. 2017;7:42229. Mu XH, Chen QW, Wu XY, et al. Gibberellins synthesis is involved in the reduction of cell flux and elemental growth rate in maize leaf under low nitrogen supply. Environ Exp Bot. 2018;150:198–208. Wenzel CL, Williamson RE, Wasteneys GO, et al. Gibberellin-induced changes in growth anisotropy precede gibberellin-dependent changes in cortical microtubule orientation in developing epidermal cells of barley leaves. Kinematic and cytological studies on a gibberellin-responsive dwarf mutant, M489. Plant Physiol. 2000;124:813–22. Nelissen H, Rymen B, Jikumaru Y, et al. A local maximum in gibberellin levels regulates maize leaf growth by spatial control of cell division. Curr Biol. 2012;22:1183–7. Zhang XJ, He LL, Zhao BL, et al. Dwarf and Increased Branching 1 controls plant height and axillary bud outgrowth in Medicago truncatula . J Exper Bot. 2020;71:6355–65. Lester DR, Ross JJ, Davies PJ, et al. Mendel’s stem length gene (Le) encodes a gibberellin 3 beta-hydroxylase. Plant Cell. 1997;9:1435–43. Itoh H, Ueguchi-Tanaka M, Sentoku N, et al. Cloning and functional analysis of two gibberellin 3β-hydroxylase genes that are differently expressed during the growth of rice. PNAS. 2001;98:8909–14. Dalmadi A, Kaló P, Jakab J, et al. Dwarf plants of diploid Medicago sativa carry a mutation in the gibberellin 3-beta-hydroxylase gene. Plant Cell Rep. 2008;27:1271–9. Hedden P, Thomas SG. Gibberellin biosynthesis and its regulation. Biochem J. 2012;444:11–25. Chen Y, Hou M, Liu L, et al. The maize DWARF1 encodes a gibberellin 3-oxidase and is dual localized to the nucleus and cytosol. Plant Physiol. 2014;166:2028–39. Phillips AL, Huttly AK, Alarcón-Reverte R, et al. GIBBERELLIN 3-OXIDASE genes regulate height and grain size in bread wheat. J Exp Bot. 2025;76:3345–58. Additional Declarations No competing interests reported. Supplementary Files 20260419SupplementaryFig.docx 20260312SupplementaryTable.xlsx Cite Share Download PDF Status: Under Review Version 1 posted Reviewers agreed at journal 15 May, 2026 Reviews received at journal 10 May, 2026 Reviewers agreed at journal 04 May, 2026 Reviewers invited by journal 03 May, 2026 Editor assigned by journal 03 May, 2026 Editor invited by journal 29 Apr, 2026 Submission checks completed at journal 29 Apr, 2026 First submitted to journal 29 Apr, 2026 You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-9460163","acceptedTermsAndConditions":true,"allowDirectSubmit":false,"archivedVersions":[],"articleType":"Research Article","associatedPublications":[],"authors":[{"id":637199722,"identity":"837c2a53-3b48-446f-936e-a7482a5f6cd7","order_by":0,"name":"Jianhang Zhang","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAABDklEQVRIie3RsWqEMBjA8U8EXT5wjZReXyFFsFc8uFdJKOhyHHRroUNEyC1CX6WPUAk45XbhFp269ja7lKqFUkoExw75g1niD/MZAJvtH7Y95HnLeoKe61btrw13llBUBe3KzSrw5R0FtoSQTIatl0ZhqWOyiNwALyhDxV8aFj9+9AoC/1gTeEi48I+vJnIrqmGW9UTSEzIFYblPCeiMC9wzE4Eq//lKfYKB0AZj4kjFBUFqJMqRhHkj4fK+H8h2Ip/zhNYTGcdXHowHo2QkYp6EpVNQ/v2T3QtMMyR6F61ZnUUSd0YSXL11XT9dZXA+95tkFRz0dfP+lFw++9o8y99wWtnweIvet9lsNpupLyy6YoZI5vLcAAAAAElFTkSuQmCC","orcid":"","institution":"Shaoxing University","correspondingAuthor":true,"prefix":"","firstName":"Jianhang","middleName":"","lastName":"Zhang","suffix":""},{"id":637199725,"identity":"c79a0345-e015-4d21-8a87-0da761f3c441","order_by":1,"name":"jiecheng Li","email":"","orcid":"","institution":"East China Normal University","correspondingAuthor":false,"prefix":"","firstName":"jiecheng","middleName":"","lastName":"Li","suffix":""},{"id":637199728,"identity":"a9a74a46-bd49-450a-b4ba-4b487526bfd8","order_by":2,"name":"Hai Xing","email":"","orcid":"","institution":"Shaoxing Cash Crop Technology Extension Center","correspondingAuthor":false,"prefix":"","firstName":"Hai","middleName":"","lastName":"Xing","suffix":""},{"id":637199731,"identity":"e6fcb938-f409-4f90-b1a2-a53e1ad95754","order_by":3,"name":"Feng Zhang","email":"","orcid":"","institution":"Nanchang University","correspondingAuthor":false,"prefix":"","firstName":"Feng","middleName":"","lastName":"Zhang","suffix":""},{"id":637199735,"identity":"dc1f8248-6de1-4ab6-ac47-ac58e10a5ca1","order_by":4,"name":"Chao Wang","email":"","orcid":"","institution":"Shaoxing University","correspondingAuthor":false,"prefix":"","firstName":"Chao","middleName":"","lastName":"Wang","suffix":""},{"id":637199740,"identity":"6adbfd75-1041-4cca-9747-de7fc3068dd4","order_by":5,"name":"Fangping Tang","email":"","orcid":"","institution":"Shaoxing University","correspondingAuthor":false,"prefix":"","firstName":"Fangping","middleName":"","lastName":"Tang","suffix":""},{"id":637199744,"identity":"0d707dee-0c21-417c-b06d-598dc7ea201e","order_by":6,"name":"Hongqing Li","email":"","orcid":"","institution":"East China Normal University","correspondingAuthor":false,"prefix":"","firstName":"Hongqing","middleName":"","lastName":"Li","suffix":""},{"id":637199745,"identity":"94f8a3df-1639-479d-b972-0ccf9e6461e3","order_by":7,"name":"Shuai Liao","email":"","orcid":"","institution":"Key Laboratory of National Forestry and Grassland Administration on East China Plant Conservation and Utilization, Shanghai Chenshan Botanical Garden","correspondingAuthor":false,"prefix":"","firstName":"Shuai","middleName":"","lastName":"Liao","suffix":""}],"badges":[],"createdAt":"2026-04-19 07:39:25","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-9460163/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-9460163/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":109022826,"identity":"f0e63ceb-2cb6-4dd2-9357-e381f76f4247","added_by":"auto","created_at":"2026-05-11 19:41:14","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":206286,"visible":true,"origin":"","legend":"\u003cp\u003eStatus of sampled leaves of \u003cem\u003eS. stolonifera\u003c/em\u003e.\u003cstrong\u003e G1\u003c/strong\u003e, first-stage adult leaf of the green-leaf type; \u003cstrong\u003eG2\u003c/strong\u003e, fully adult leaf of the green-leaf type;\u003cstrong\u003e YG\u003c/strong\u003e, young leaf of the green-leaf type;\u003cstrong\u003e P1\u003c/strong\u003e, first-stage adult leaf of the purple-variegated type; \u003cstrong\u003eP2\u003c/strong\u003e, fully adult leaf of the purple-variegated type; \u003cstrong\u003eYP\u003c/strong\u003e, young leaf of the purple-variegated type. \u003cstrong\u003eW1\u003c/strong\u003e, first-stage adult leaf of the white-variegated type; \u003cstrong\u003eW2\u003c/strong\u003e, fully adult leaf of the white-variegated type;\u003cstrong\u003e YW\u003c/strong\u003e, young leaf of the white-variegated type.\u003c/p\u003e","description":"","filename":"floatimage1.png","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/51e95a91c3e198c2fb7e5845.png"},{"id":109022828,"identity":"9b809b56-6518-46b0-969a-a62cb5fc65a9","added_by":"auto","created_at":"2026-05-11 19:41:14","extension":"jpeg","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":248198,"visible":true,"origin":"","legend":"\u003cp\u003eOverview of \u003cem\u003eS. stolonifera\u003c/em\u003egenome assembly. \u003cstrong\u003e(a)\u003c/strong\u003e Photographs of \u003cem\u003eS. stolonifera\u003c/em\u003e. \u003cstrong\u003e(b)\u003c/strong\u003e Chromosome counts in \u003cem\u003eS. stolonifera\u003c/em\u003e root tip cells. \u003cstrong\u003e(c)\u003c/strong\u003e Hi-C interaction heatmap of \u003cem\u003eS. stolonifera\u003c/em\u003e. \u003cstrong\u003e(d)\u003c/strong\u003e Circos plot of the \u003cem\u003eS. stolonifera\u003c/em\u003e genome assembly. \u003cstrong\u003eⅰ\u003c/strong\u003e, syntenic blocks; \u003cstrong\u003eⅱ\u003c/strong\u003e, gene density and unknown base (N) ratio; \u003cstrong\u003eⅲ\u003c/strong\u003e, GC contentand skew; G1, W1, and P1, sequencing coverage for three phenotypes (G, green-leaf type; W, white-variegated type; P, purple-variegated type).\u003c/p\u003e","description":"","filename":"floatimage2.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/c3be792bfa3f29df92a04a9f.jpeg"},{"id":109022864,"identity":"2ede4b26-3063-4ad2-be46-a07d0257d7ee","added_by":"auto","created_at":"2026-05-11 19:41:15","extension":"jpeg","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":440532,"visible":true,"origin":"","legend":"\u003cp\u003eComparative genomic and evolutionary analysis of \u003cem\u003eS. stolonifera\u003c/em\u003e. \u003cstrong\u003e(a)\u003c/strong\u003e Copy number distribution of all gene families in 25 selected angiosperm species. \u003cstrong\u003e(b)\u003c/strong\u003e Gene family cluster petal diagram, with the central circle representing common gene families and the outer petals depicting species-specific gene families.\u003cstrong\u003e(c)\u003c/strong\u003e Phylogenetic tree of 25 plant species. Numbers represent divergence times at each node (Mya, million years ago). Gene family expansions and contractions are indicated by green and red numbers, respectively.\u003cstrong\u003e (d)\u003c/strong\u003e KS distribution analysis. \u003cstrong\u003e(e)\u003c/strong\u003e Syntenic depth analyses of \u003cem\u003eSaxifraga stolonifera \u003c/em\u003eand \u003cem\u003eVitis vinifera\u003c/em\u003e genomes.\u003c/p\u003e","description":"","filename":"floatimage3.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/9356ce8cfae49b90b8a50497.jpeg"},{"id":109068081,"identity":"ba5bf67f-b5e7-42da-beff-c2dfdf994007","added_by":"auto","created_at":"2026-05-12 10:03:28","extension":"jpeg","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":895710,"visible":true,"origin":"","legend":"\u003cp\u003eTransverse sections of different \u003cem\u003eS. stolonifera\u003c/em\u003eleaves under cold-field emission scanning electron microscopy, showing air spaces. \u003cstrong\u003e(a–c)\u003c/strong\u003e sampling position (from left to right, the green-leaf type G, white-variegated type W, purple-variegated type P, the same below). \u003cstrong\u003e(d–f)\u003c/strong\u003eleaf transverse section (130×). \u003cstrong\u003e(g–i)\u003c/strong\u003etransverse section of the leaf vein area (300×). \u003cstrong\u003e(j–l)\u003c/strong\u003e transverse section between the leaf vein area (300×). \u003cstrong\u003eEc\u003c/strong\u003e, epidermal cells. \u003cstrong\u003ePt\u003c/strong\u003e, palisade tissue. \u003cstrong\u003eVb\u003c/strong\u003e, vascular bundle. Arrows in h and i indicate air spaces.\u003c/p\u003e","description":"","filename":"floatimage4.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/4ebd2aa82cafefe947547af7.jpeg"},{"id":109022827,"identity":"b891f15f-c444-4039-8e0b-abf2600585fb","added_by":"auto","created_at":"2026-05-11 19:41:14","extension":"jpeg","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":220538,"visible":true,"origin":"","legend":"\u003cp\u003eMetabolomic profiling of anthocyanins and plant hormones.\u003cstrong\u003e (a)\u003c/strong\u003e Heatmaps of differentially abundant metabolites (DAMs) in the three compared phenotypes. \u003cstrong\u003e(b)\u003c/strong\u003e Heatmaps of plant hormones in three leaf phenotypes. \u003cstrong\u003e(c) \u003c/strong\u003eVenn analysis of differential plant hormones in four comparison groups. \u003cstrong\u003eg\u003c/strong\u003e, second-stage adult leaf of the green-leaf type; \u003cstrong\u003ep\u003c/strong\u003e, second-stage adult leaf of the purple-variegated type; \u003cstrong\u003ew\u003c/strong\u003e, second-stage adult leaf of the white-variegated type.\u003c/p\u003e","description":"","filename":"floatimage5.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/be42bdf597567f8a2567618a.jpeg"},{"id":109022825,"identity":"91201d62-9373-4588-8313-ef0d5a26cf06","added_by":"auto","created_at":"2026-05-11 19:41:14","extension":"jpeg","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":534067,"visible":true,"origin":"","legend":"\u003cp\u003eAnthocyanin biosynthetic pathways and gene expression in three leaf phenotypes. Gene expression levels (log\u003csub\u003e10\u003c/sub\u003e (TPM + 1)) in the three leaf phenotypes are represented by color grading, where TPM = the average value of each group of all samples. \u003cstrong\u003eG\u003c/strong\u003e, green-leaf type; \u003cstrong\u003eP\u003c/strong\u003e, purple-variegated type; \u003cstrong\u003eW\u003c/strong\u003e, white-variegated type.\u003c/p\u003e","description":"","filename":"floatimage6.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/bb02b1a5458de329d9156094.jpeg"},{"id":109022862,"identity":"f98a572f-b633-4794-936b-b1f2e9286a5e","added_by":"auto","created_at":"2026-05-11 19:41:15","extension":"jpeg","order_by":7,"title":"Figure 7","display":"","copyAsset":false,"role":"figure","size":461069,"visible":true,"origin":"","legend":"\u003cp\u003eTranscriptomic and plant hormonecorrelation analysis in three leaf phenotypes. \u003cstrong\u003e(a)\u003c/strong\u003e Map and heatmaps of the gibberellin biosynthesis pathway. \u003cstrong\u003e(b) \u003c/strong\u003eDendrogram showing co-expression modules identified by WGCNA across plant hormones in three leaf phenotypes. \u003cstrong\u003e(c)\u003c/strong\u003eHeatmap displaying the correlations between gene expression modules and gibberellin biosynthesis. Each of the 23 rows corresponds to a specific module (ME) indicated by a distinct color.\u003c/p\u003e","description":"","filename":"floatimage7.jpeg","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/817c61e0425796e097f2ef7c.jpeg"},{"id":109069329,"identity":"42937563-154a-4bcf-9cb0-53df12858989","added_by":"auto","created_at":"2026-05-12 10:22:37","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":3437702,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/50286f29-722e-4009-beb8-35d4a9ec2cea.pdf"},{"id":109022859,"identity":"1c331156-998e-4072-ac75-f99c2040f28c","added_by":"auto","created_at":"2026-05-11 19:41:14","extension":"docx","order_by":0,"title":"","display":"","copyAsset":false,"role":"supplement","size":1248569,"visible":true,"origin":"","legend":"","description":"","filename":"20260419SupplementaryFig.docx","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/34a51ae53f54f08913a5caf5.docx"},{"id":109022824,"identity":"a5e79640-6962-47e4-ad15-8df30f6c02a8","added_by":"auto","created_at":"2026-05-11 19:41:14","extension":"xlsx","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":8052991,"visible":true,"origin":"","legend":"","description":"","filename":"20260312SupplementaryTable.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-9460163/v1/6540784a1e64f8b0b55c709d.xlsx"}],"financialInterests":"No competing interests reported.","formattedTitle":"Multi-omics reveals divergent regulation of anthocyanin glycosylation and gibberellin biosynthesis underlying leaf variegation in Saxifraga stolonifera","fulltext":[{"header":"Introduction","content":"\u003cp\u003eVariegated leaf plants represent a distinct group of colored-leaf species characterized by stable, heritable patterns formed by differentially pigmented leaf regions. These features make them valuable model systems for studying plant chromatology, mechanisms of adaptation to low-light environments, and chloroplast development[\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e, \u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e]. These variegated leaves can be classified into five types: chlorophyll type, air space type, epidermis type, pigment type, and appendages type[\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e]. Among these, the air space type and pigment type are the most prevalent and often co-occur on the same leaf, in species such as \u003cem\u003eActinidia lanceolata\u003c/em\u003e, \u003cem\u003eErythronium dens-canis\u003c/em\u003e, and \u003cem\u003eSaxifraga stolonifera\u003c/em\u003e[\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e, \u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e]. Variegated leaves contribute to multiple adaptive functions, including acclimation to shade conditions, enhanced cold tolerance, defense against herbivores, and facilitation of self-pollination[\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e, \u003cspan additionalcitationids=\"CR6 CR7\" citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e]. In addition, variegated plants possess significant ornamental value due to their visually striking foliage.\u003c/p\u003e \u003cp\u003eAlthough naturally occurring leaf variegation has been studied for decades, previous research has largely focused on morphological classification and ecological functions; genomics, transcriptomics, and other omics studies have been limited by comparison. To date, transcriptomic studies of variegated leaves have primarily focused on the types of variegations associated with chlorophyll deficiency or red-leaf phenotypes. These investigations have elucidated the molecular mechanisms underlying chlorophyll degradation or anthocyanin accumulation in variegated tissues, as demonstrated in \u003cem\u003eAnanas comosus\u003c/em\u003e var. \u003cem\u003ebracteatus\u003c/em\u003e[\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e], \u003cem\u003eEpipremnum aureum\u003c/em\u003e \u0026lsquo;Marble Queen\u0026rsquo;[\u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e], \u003cem\u003eClivia miniata\u003c/em\u003e var. \u003cem\u003evariegata\u003c/em\u003e[\u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e], \u003cem\u003eAlternanthera bettzickiana\u003c/em\u003e[\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e], \u003cem\u003eCamellia sinensis\u003c/em\u003e \u0026lsquo;Zijuan\u0026rsquo;[\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e], and \u003cem\u003ePerilla frutescens\u003c/em\u003e[\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e]. In transcriptome analyses of air space-type variegation in \u003cem\u003eTrifolium pratense\u003c/em\u003e and \u003cem\u003ePrimulina pungentisepala\u003c/em\u003e, key pathways related to photosynthesis, redox regulation, cell-wall modification, and nitrogen metabolism have been examined[\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e, \u003cspan citationid=\"CR15\" class=\"CitationRef\"\u003e15\u003c/span\u003e]. Anthocyanins, a major class of water-soluble flavonoid pigments, play a critical role in determining the coloration of plant organs and tissues[\u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e16\u003c/span\u003e, \u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e18\u003c/span\u003e]. Under high-light conditions, the flavonoid biosynthesis pathway is activated, resulting in the accumulation of photoprotective and antioxidative flavonoids, particularly flavonols and anthocyanins, in leaf tissues[\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e]. The spatially restricted accumulation of anthocyanins is primarily responsible for pigment-type variegation[\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e, \u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e20\u003c/span\u003e]. For example, in \u003cem\u003eCorydalis hemidicentra\u003c/em\u003e, insertion of a 254-bp transposon into the \u003cem\u003ebHLH35\u003c/em\u003e gene enhances anthocyanin biosynthesis, conferring an environmentally adaptive gray phenotype[\u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e]. Despite these advances, the developmental dynamics and transcriptional regulatory networks governing both air space-type and pigment-type variegation remain poorly understood[\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e].\u003c/p\u003e \u003cp\u003e \u003cem\u003eSaxifraga stolonifera\u003c/em\u003e is a shade-tolerant perennial exhibiting variegated leaves with both ornamental and medicinal properties, and it displays diverse leaf variegation patterns that support adaptation across variable environments[\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e, \u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e, \u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e]. In its natural habitat, this species commonly shows three stable leaf phenotypes: green-leaf type, white-variegated type, and purple-variegated type. Here, we report a high-quality genome assembly of \u003cem\u003eS. stolonifera\u003c/em\u003e generated by integrating Illumina short-read sequencing, PacBio long-read sequencing, and high-throughput chromosome conformation capture (Hi-C) technologies. Furthermore, we performed comprehensive metabolomics and transcriptomic profiling of young leaves, first-stage adult leaves, and fully mature leaves across the three phenotypes. A weighted gene co-expression network analysis (WGCNA) was also conducted to identify key genes and regulatory modules associated with changes in hormone-mediated metabolites during leaf variegation. By integrating morphological and physiological traits with plant hormone profiles, metabolomics data, and transcriptomic landscapes, this study provides a high-quality reference genome and offers in-depth insights into the molecular mechanisms underlying leaf variegation. These resources establish a foundation for future genetic, genomic, and functional studies within the genus \u003cem\u003eSaxifraga\u003c/em\u003e.\u003c/p\u003e"},{"header":"Materials and methods","content":"\u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003ePlant materials\u003c/h2\u003e \u003cp\u003eThe green-leaf (G) and purple-variegated (P) phenotypes of \u003cem\u003eS. stolonifera\u003c/em\u003e were collected from Tianmu Mountain, Zhejiang, China (30\u0026deg;18'30\" N, 119\u0026deg;23'47\u0026Prime; E), while the white-variegated (W) phenotype was obtained from Doupeng Mountain, Guizhou, China (26\u0026deg;20\u0026prime;08\u0026Prime; N, 107\u0026deg;17\u0026prime;34\u0026Prime; E). The formal identification of \u003cem\u003eS. stolonifera\u003c/em\u003e was performed by Dr. Jianhang Zhang. Voucher specimens were deposited in herbarium of East China Normal University (HSNU). As S. \u003cem\u003estolonifera\u003c/em\u003e is not a protected species, no collection permit was required. Starting in September 2020, all three phenotypes were cultivated in a greenhouse at the Biological Station of East China Normal University using a 1:3 mixture of nutrient soil and vermiculite as substrate. On April 20, 2021, healthy, uniformly growing individuals of each phenotype were selected for experimentation. From each plant, two leaves were sampled: one first-stage adult leaf and one fully developed adult leaf of similar size and growth status. The smallest leaf exhibiting typical morphological characteristics was designated as the first-stage adult leaf. On June, 22, 2025, the young leaves lacking characteristic pigmentation patterns were sampled. In this study, the second emerging leaf was defined as the young leaf, the third as the first-stage adult leaf (showing distinct patterning), and the sixth or seventh leaf as the mature adult leaf. Leaf developmental stages are illustrated in Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003e. After harvest, leaves were wrapped in aluminum foil, rapidly frozen in liquid nitrogen, and stored at \u0026minus;\u0026thinsp;80\u0026deg;C for downstream analyses. In April 2023, individuals of all three phenotypes were transplanted from East China Normal University to the greenhouse at Shaoxing University. On May, 7, 2024, root tips were collected for chromosome karyotype analysis[\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e], and whole plants of the green-leaf phenotype at the flowering stage were sampled for genome sequencing.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e\n\u003ch3\u003eMacroscopic leaf morphology and phenotypic stability\u003c/h3\u003e\n\u003cp\u003eEighteen plants, six representing each of the three leaf phenotypes with uniform growth characteristics, were selected for morphological analysis. To assess the adaptability of the green-leaf and purple-variegated types to varying light conditions, these individuals were exposed to distinct light intensities within the same greenhouse: full sunlight, shading with a single-layer shade net (reducing light intensity to approximately 25% of full sunlight), and shading with a double-layer shade net (reducing light intensity approximately 10%). Substrate composition, watering regimen frequency, and fertilization protocols were standardized across all experimental groups. Weekly visual documentation of leaf morphology and color variation was performed using a Canon EOS 800D camera, while leaf surface features were examined using a KEYENCE VHX-5000 stereoscope (Japan).\u003c/p\u003e\n\u003ch3\u003eLeaf structure and ultrastructure\u003c/h3\u003e\n\u003cp\u003eFully developed leaves of the three \u003cem\u003eS. stolonifera\u003c/em\u003e phenotypes exhibiting healthy and uniform growth were used for structural and ultrastructural analyses. Leaves were manually sectioned, and thin, intact slices were carefully selected to prepare temporary water-sealed slides. Cell morphology, dimensions, and arrangement; chloroplast number, distribution, and spatial organization; and the spatial arrangement of anthocyanin-containing cells were observed under bright-field mode at 4 \u0026times;, 10 \u0026times;, and 20 \u0026times; magnification using an Echo RVL-100-M inverted integrated fluorescence microscope (Discover ECHO, US). Chloroplast autofluorescence was assessed in fluorescence mode via the CY5 channel to determine chloroplast abundance and intracellular localization. For scanning electron microscopy, leaf blades were longitudinally dissected along both sides of the midvein, preserving a central segment approximately 1-2cm wide to include both regular and variegated regions. Observations were conducted using a Hitachi S-4800 cold-field emission scanning electron microscope (Japan), operated at 1 kV.\u003c/p\u003e\n\u003ch3\u003eMeasurements of the content of chlorophyll and anthocyanins\u003c/h3\u003e\n\u003cp\u003eThe chlorophyll content was measured using a SPAD-502 Plus meter (Konica Minolta, Japan) and expressed as relative SPAD values. A total of 480 leaves were sampled from 40 uniformly growing plants representing the three leaf phenotypes. One leaf each from the first-stage adult and fully adult stages, matched for size and developmental status, was selected per plant. For fully adult leaves, measurements were taken either over the main vein (white-variegated area) or interveinal areas (purple-variegated area). Due to the difficulty of distinguishing variegated and non-variegated areas on first-stage adult leaves, measurements were restricted to the main vein area. Data were processed using Microsoft Excel (2019), and differences in the chlorophyll content among leaf phenotypes were analyzed by one-way ANOVA using R v4.1.2.\u003c/p\u003e \u003cp\u003eThe anthocyanin content was quantified using the MetWare platform (Wuhan, China; \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://www.metware.cn/\u003c/span\u003e\u003cspan address=\"http://www.metware.cn/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e) and the AB Sciex QTRAP 6500 LC-MS/MS system. Nine \u003cem\u003eS. stolonifera\u003c/em\u003e plants, three per phenotype, with similar growth were selected. One first-stage adult leaf and one fully adult leaf, matched for size and developmental stage, were collected from each plant, yielding nine samples in total. Samples were immediately frozen in liquid nitrogen and stored at -80 ℃ until analysis. The results of hierarchical clustering analysis (HCA) of both samples and metabolites were visualized as dendrograms accompanied by heatmaps. HCA was performed using the R package pheatmap with normalized metabolite signal intensities (unit variance scaling) represented on a color scale. Differentially abundant metabolites (DAMs) were identified based on pairwise comparisons between groups: adult leaves of the purple-variegated type (p) versus green-leaf type (g); white-variegated type (w) versus g; and p versus w. DAMs were considered significant if they met the following criteria: absolute Log\u003csub\u003e2\u003c/sub\u003e(FC)\u0026thinsp;\u0026ge;\u0026thinsp;1 (i.e., fold change\u0026thinsp;\u0026ge;\u0026thinsp;2 or\u0026le; -2), P-value\u0026thinsp;\u0026le;\u0026thinsp;0.05, and presence rate\u0026thinsp;\u0026ge;\u0026thinsp;1 across all samples. Identified metabolites were annotated using the KEGG COMPOUND database (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://www.kegg.jp/kegg/compound/\u003c/span\u003e\u003cspan address=\"http://www.kegg.jp/kegg/compound/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e).\u003c/p\u003e\n\u003ch3\u003eDetection of the phytohormone content\u003c/h3\u003e\n\u003cp\u003ePhytohormone levels were quantified using an ultra-high-performance liquid chromatography-tandem mass spectrometry (UHPLC-MS/MS) system (ExionLC\u0026trade; AD UHPLC-QTRAP 6500+, AB SCIEX Corp., Boston, MA, USA) at Novogene Co., Ltd. (Beijing, China). On 22 June 2025, young leaves and first-stage adult leaves were sampled. Nine healthy plants, three per phenotype (G, P, W), with consistent growth were selected. From each plant, one young leaf (second leaf) and one first-stage adult leaf (third leaf), comparable in size and developmental stage, were harvested, resulting in 18 samples. All samples were flash-frozen in liquid nitrogen and stored at -80℃ prior before analysis. HCA was performed on all detected plant hormones, with hormone concentrations normalized and clustered accordingly. Differentially abundant metabolites were identified through pairwise comparisons using the thresholds of fold change (FC)\u0026thinsp;\u0026gt;\u0026thinsp;1.2 or FC\u0026thinsp;\u0026lt;\u0026thinsp;0.833 and P-value\u0026thinsp;\u0026lt;\u0026thinsp;0.05.\u003c/p\u003e \u003cdiv id=\"Sec8\" class=\"Section2\"\u003e \u003ch2\u003eGenome assembly and annotation\u003c/h2\u003e \u003cp\u003eGenome sequencing was conducted using the Illumina platform at Biomarker Technologies (Beijing, China), integrating PacBio HiFi, Hi-C, and RNA-seq data. For genome annotation, mixed RNA samples were collected from four healthy tissues (root, stem, leaf, and flower) and subjected to RNA sequencing. Genomic DNA was extracted from \u003cem\u003eS. stolonifera\u003c/em\u003e leaves to construct 150-bp paired-end libraries, which were sequenced on the Illumina HiSeq X platform. Circular consensus sequencing (CCS) reads generated by PacBio HiFi were assembled using Hifiasm (version 0.19.9-r616; \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://github.com/chhylp123/hifiasm\u003c/span\u003e\u003cspan address=\"https://github.com/chhylp123/hifiasm\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e) with key parameters \u0026lsquo;-n 3 -I\u0026thinsp;=\u0026thinsp;0\u0026rsquo; and default settings. Chromosome-level scaffolding was performed using HapHiC (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://github.com/zengxiaofei/HapHiC\u003c/span\u003e\u003cspan address=\"https://github.com/zengxiaofei/HapHiC\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e) with default parameters. Genome completeness was evaluated using BUSCO (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://github.com/metashot/busco\u003c/span\u003e\u003cspan address=\"https://github.com/metashot/busco\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e) in genome mode (-m genome). Jellyfish v2.3.05[\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e] was used to count k-mers (K\u0026thinsp;=\u0026thinsp;31), and GCE v1.0.2[\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e] was applied to estimate genome size and heterozygosity from k-mer frequency distributions.\u003c/p\u003e \u003cp\u003eTransposable elements (TEs) were identified using a combined homology-based and de novo approach. A de novo repeat library was constructed using RepeatModeler2 v2.0.1[\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e], which incorporates RECON v1.0.8[\u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e] and RepeatScout v1.0.6[\u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e]. Full-length LTR retrotransposons were initially detected using LTRharvest v1.5.10[\u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e30\u003c/span\u003e] and LTR_finder v1.07[\u003cspan citationid=\"CR31\" class=\"CitationRef\"\u003e31\u003c/span\u003e], and then refined with LTR_retriever. A non-redundant, species-specific TE library was generated by merging the de novo library with the Dfam v3.5 database. Final TE sequences were classified via RepeatMasker v4.12[\u003cspan citationid=\"CR32\" class=\"CitationRef\"\u003e32\u003c/span\u003e] using homology searches. Tandem repeats were annotated using TRF v4.09[\u003cspan citationid=\"CR33\" class=\"CitationRef\"\u003e33\u003c/span\u003e] and MISA v2.1[\u003cspan citationid=\"CR34\" class=\"CitationRef\"\u003e34\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eGene annotation was performed using the Makerpipeline[\u003cspan citationid=\"CR35\" class=\"CitationRef\"\u003e35\u003c/span\u003e]. To minimize annotation errors, gene predictions were performed using three complementary approaches: de novo prediction, transcript-based assembly, and homology-based alignment. Homology-based predictions utilized protein sequences from \u003cem\u003eArabidopsis thaliana\u003c/em\u003e, \u003cem\u003eBergenia scopulosa\u003c/em\u003e, \u003cem\u003eChrysosplenium sinicum\u003c/em\u003e, \u003cem\u003eLiquidambar formosana\u003c/em\u003e, and \u003cem\u003eVitis vinifera\u003c/em\u003e, which were mapped to the \u003cem\u003eS. stolonifera\u003c/em\u003egenome using GeMoMa v1.7[\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e]. De novo gene prediction was carried out using Augustus v3.2.3[\u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e] and SNAP[\u003cspan citationid=\"CR38\" class=\"CitationRef\"\u003e38\u003c/span\u003e]. RNA-seq reads were aligned and assembled using HISAT2 v2.1.0[\u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e] and StringTie v2.1.4[\u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e40\u003c/span\u003e], respectively, and used as evidence for gene prediction with GeneMarkS-T v5.1[\u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e41\u003c/span\u003e]. Additionally, PASA v2.4.1[\u003cspan citationid=\"CR42\" class=\"CitationRef\"\u003e42\u003c/span\u003e] was employed to predict gene models from RNA-Bloom v2.0.0[\u003cspan citationid=\"CR43\" class=\"CitationRef\"\u003e43\u003c/span\u003e] assembled unigenes and PacBio/ONT full-length transcripts. All predicted models were consolidated using EVM v1.1.1[\u003cspan citationid=\"CR44\" class=\"CitationRef\"\u003e44\u003c/span\u003e] and subsequently updated with PASA. The completeness of the annotated gene set was evaluated using BUSCO v5.2.2[\u003cspan citationid=\"CR45\" class=\"CitationRef\"\u003e45\u003c/span\u003e] in protein mode (-m proteins).\u003c/p\u003e \u003c/div\u003e\n\u003ch3\u003eComparative genomic analysis\u003c/h3\u003e\n\u003cp\u003eA comparative genomic analysis was conducted involving \u003cem\u003eS. stolonifera\u003c/em\u003e and 24 other plant species (Table \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e). Protein sequences were clustered into gene families using OrthoFinder v2.4.0[\u003cspan citationid=\"CR46\" class=\"CitationRef\"\u003e46\u003c/span\u003e]. Functional annotations of gene families were derived from the Pfam V33.1 database[\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e]. Unique gene families in each species were identified through GO and KEGG enrichment analyses. Single-copy orthologous genes present in at least 80.0% of species (n\u0026thinsp;=\u0026thinsp;799) were aligned using MAFFT v7.205 (--localpair --maxiterate 1000)[\u003cspan citationid=\"CR48\" class=\"CitationRef\"\u003e48\u003c/span\u003e]. A maximum-likelihood phylogenetic tree was inferred using iqtree v2.2.0 (JTT\u0026thinsp;+\u0026thinsp;F+I+G4, 1000 bootstrap replicates)[\u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e49\u003c/span\u003e], with \u003cem\u003eAmborella trichopoda\u003c/em\u003e designated as the outgroup. Divergence times were estimated using the MCMCTree program in the PAML v4.9i package under default settings. Gene family expansions and contractions (family-wide P-values and Viterbi P-values\u0026thinsp;\u0026lt;\u0026thinsp;0.05) were analyzed using CAF\u0026Eacute; v4.2.1[\u003cspan citationid=\"CR50\" class=\"CitationRef\"\u003e50\u003c/span\u003e]. The evolutionary dynamics of the \u003cem\u003eS. stolonifera\u003c/em\u003e genome were further explored by calculating the synonymous substitution rate (Ks) for collinear gene pairs using WGDI v0.71[\u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e51\u003c/span\u003e].\u003c/p\u003e\n\u003ch3\u003eRNA sequencing and bioinformatic analysis\u003c/h3\u003e\n\u003cp\u003eTotal RNA was extracted from 38 samples, including one blind test sample from each of the P1 and W2 groups, comprising first-stage adult and fully adult leaves of the three leaf phenotypes. RNA purity and integrity were assessed using a NanoDrop ND-1000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA) and the Agilent RNA Nano 6000 Assay Kit on the Agilent Bioanalyzer 2100 (Agilent Technologies, Santa Clara, CA, USA). Qualified libraries were subjected to next-generation sequencing (NGS) on the Illumina platform at Biomarker Technologies (Beijing, China). Additional transcriptome data for young leaves, representing the three phenotypes before characteristic patterns become visible, were collected in 2025 and sequenced on the T7 platform at Novogene Co., Ltd. (Beijing, China).\u003c/p\u003e \u003cp\u003eClean reads were aligned to the \u003cem\u003eS. stolonifera\u003c/em\u003e reference genome using HISAT2[\u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e]. Read counts were quantified using the subread package in R (featureCounts), and gene expression levels were estimated using TPM (transcripts per million) via StringTie[\u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e40\u003c/span\u003e]. In the KEGG database, KO (KEGG Orthology) identifiers denote functionally conserved orthologous gene groups. Based on whole-genome functional annotation, all genes associated with the phenylpropanoid biosynthesis (ko00940), flavonoid biosynthesis (ko00941), and anthocyanin biosynthesis (ko00942) pathways were identified[\u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e16\u003c/span\u003e, \u003cspan additionalcitationids=\"CR53 CR54\" citationid=\"CR52\" class=\"CitationRef\"\u003e52\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR55\" class=\"CitationRef\"\u003e55\u003c/span\u003e]. Genes involved in diterpenoid biosynthesis (map00904) were also annotated [\u003cspan citationid=\"CR56\" class=\"CitationRef\"\u003e56\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eProtein domains of candidate genes were identified using HMMER via the Quick Gene Family Identification plugin in TBtools-II[\u003cspan citationid=\"CR57\" class=\"CitationRef\"\u003e57\u003c/span\u003e], and genes with inconsistent domain architectures were filtered. DEGs were identified from transcriptome data using DESeq2[\u003cspan citationid=\"CR58\" class=\"CitationRef\"\u003e58\u003c/span\u003e] with the following criteria: false discovery rate (FDR)\u0026thinsp;\u0026lt;\u0026thinsp;0.05 and |log\u003csub\u003e2\u003c/sub\u003e fold change| \u0026ge; 1. The Venny tool (accessed 27 October 2025; \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.omicshare.com/tools/\u003c/span\u003e\u003cspan address=\"https://www.omicshare.com/tools/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e) was used to visualize overlapping DEGs. After discarding genes with low relative expression (FPKM\u0026thinsp;\u0026gt;\u0026thinsp;1 in more than 90% of the samples), the WGCNA plugin (Format =normalized count, Normalized method\u0026thinsp;=\u0026thinsp;raw, Sample percentage\u0026thinsp;=\u0026thinsp;0.9, Expression Cutoff\u0026thinsp;=\u0026thinsp;1, Filter Method\u0026thinsp;=\u0026thinsp;MAD, Reserved genes Num. =20,000; R\u003csup\u003e2\u003c/sup\u003e coutoff\u0026thinsp;=\u0026thinsp;0.8, Recommended\u0026thinsp;=\u0026thinsp;6, Customized\u0026thinsp;=\u0026thinsp;8; min Module Size\u0026thinsp;=\u0026thinsp;30, module cuttree height\u0026thinsp;=\u0026thinsp;0.25, select max blocksize\u0026thinsp;=\u0026thinsp;20,000; x axis label angle\u0026thinsp;=\u0026thinsp;45, KME cutoff\u0026thinsp;=\u0026thinsp;0.2, Choose method\u0026thinsp;=\u0026thinsp;2) in TBtools-II was used to identify modules of highly correlated genes attributable to plant hormones based on the filtered FPKM data[\u003cspan citationid=\"CR57\" class=\"CitationRef\"\u003e57\u003c/span\u003e, \u003cspan citationid=\"CR59\" class=\"CitationRef\"\u003e59\u003c/span\u003e].\u003c/p\u003e"},{"header":"Results","content":"\u003cp\u003e \u003cb\u003eGenome assembly and annotation of\u003c/b\u003e \u003cb\u003eS. stolonifera\u003c/b\u003e\u003c/p\u003e \u003cp\u003eTo establish a genomic framework for dissecting leaf variegation mechanisms, we generated a high-quality chromosome-scale genome assembly of S. stolonifera. A green-phenotype individual of \u003cem\u003eS. stolonifera\u003c/em\u003e from the Yonglai village population on Qingliangfeng Mountain in Jixi County, Anhui Province, China, was selected for sequencing (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ea). Karyotype analysis confirmed diploidy with 2n\u0026thinsp;=\u0026thinsp;36 chromosomes (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eb). K-mer analysis estimated a genome size of ~\u0026thinsp;1.63 Gb with 0.29% heterozygosity and 57.6% repetitive sequence, based on the 31-mer frequency distribution (Fig. \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e). To achieve a high-quality genome assembly, we integrated multiple sequencing technologies by integrating Illumina short-read (116.53 Gb), PacBio HiFi long-read (99.65 Gb), and Hi-C (325.13 Gb) sequencing data (Table \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e). The initial assembly yielded 2,139 contigs, which were subsequently anchored into 18 pseudo-chromosomes using Hi-C interaction maps, accounting for 94% of the assembly (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ec, d; Tables\u0026nbsp;\u003cspan refid=\"Tab1\" class=\"InternalRef\"\u003e1\u003c/span\u003e, S1, S2). The final genome assembly spans 2.01 Gb with a scaffold N50 of 3,535,068 bp and a BUSCO completeness score of 98.5%, indicating high contiguity and completeness (Table\u0026nbsp;\u003cspan refid=\"Tab1\" class=\"InternalRef\"\u003e1\u003c/span\u003e; Table S3).\u003c/p\u003e \u003cp\u003e \u003cdiv class=\"gridtable\"\u003e\u003ctable float=\"Yes\" id=\"Tab1\" border=\"1\"\u003e \u003ccaption language=\"En\"\u003e \u003cdiv class=\"CaptionNumber\"\u003eTable 1\u003c/div\u003e \u003cdiv class=\"CaptionContent\"\u003e \u003cp\u003eStatistics for the final genome assembly of \u003cem\u003eS. stolonifera\u003c/em\u003e\u003c/p\u003e \u003c/div\u003e \u003c/caption\u003e \u003ccolgroup cols=\"2\"\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c1\" colnum=\"1\"\u003e\u003c/div\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c2\" colnum=\"2\"\u003e\u003c/div\u003e \u003cthead\u003e \u003ctr\u003e \u003cth align=\"left\" colname=\"c1\"\u003e \u003cp\u003eGenome information\u003c/p\u003e \u003c/th\u003e \u003cth align=\"left\" colname=\"c2\"\u003e \u003cp\u003e\u003cem\u003eS. stolonifera\u003c/em\u003e (PacBio\u0026thinsp;+\u0026thinsp;Hi-C)\u003c/p\u003e \u003c/th\u003e \u003c/tr\u003e \u003c/thead\u003e \u003ctbody\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003ePacBio HiFi data (Gb)\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e99.66\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003eHi-C clean data (Gb)\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e325.13\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003eAssembly size (Gb)\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e2.01\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003eNumber of contigs\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e2,139\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003eScaffold N50 (bp)\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e3,535,068\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003eScaffold max (bp)\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e14,862,665\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003eAnchor ratio (%)\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e94.00%\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003eGC content\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e34.12%\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cp\u003eBUSCO (%)\u003c/p\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cp\u003e98.50%\u003c/p\u003e \u003c/td\u003e \u003c/tr\u003e \u003c/tbody\u003e \u003c/colgroup\u003e \u003c/table\u003e\u003c/div\u003e \u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eRepetitive sequences constitute a major fraction of eukaryotic genomes and primarily include tandem repeats and interspersed repeats[\u003cspan citationid=\"CR60\" class=\"CitationRef\"\u003e60\u003c/span\u003e, \u003cspan citationid=\"CR61\" class=\"CitationRef\"\u003e61\u003c/span\u003e]. To annotate repetitive elements in the \u003cem\u003eS. stolonifera\u003c/em\u003e genome, both de novo and homology-based prediction methods were employed. A total of 1,835,418 repetitive sequences were identified, spanning 822,048,413 bp (~\u0026thinsp;73.49% of the assembled genome) (Table S4). Retroelements were the most abundant repeat class, accounting for 44.99% of the genome, followed by DNA transposons (14%) and tandem repeats (14.5%) (Table S4).\u003c/p\u003e \u003cp\u003eGene models for the \u003cem\u003eS. stolonifera\u003c/em\u003e genome were predicted through an integrative approach combining ab initio, homology-based, and transcriptome-supported predictions. This strategy yielded a total of 37,191 protein-coding genes, with an average gene length of 3,914.91 bp and an average of 4.95 exons per gene (Fig. \u003cspan refid=\"MOESM2\" class=\"InternalRef\"\u003eS2\u003c/span\u003e; Table S5). Of these, 36,664 genes were successfully annotated, resulting in an annotation rate of 98.58% (Tables S6, S7). TrEMBL analysis indicated that 36,585 genes (98.37%) were functionally annotated, with 28,614 encoding metabolic enzymes as classified by Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways (Table S7). The completeness and quality of the annotation were evaluated using BUSCO. Among the 1,614 conserved orthologs assessed, 1588 (98.39%) were complete, including 86.18% single-copy and 12.21% duplicated genes, indicating high annotation completeness (Table S8). This high-quality genome provides a robust reference for subsequent functional and evolutionary analyses.\u003c/p\u003e \u003cp\u003e \u003cb\u003eGenome evolutionary and whole-genome duplication analysis of\u003c/b\u003e \u003cb\u003eS. stolonifera\u003c/b\u003e\u003c/p\u003e \u003cp\u003eTo place \u003cem\u003eS. stolonifera\u003c/em\u003e in an evolutionary context, we compared its genome with 24 representative plant species (Table S9). A total of 4,130 shared gene families were identified across all species, while 402 gene families were unique to \u003cem\u003eS. stolonifera\u003c/em\u003e, enriched in pathways related to aminoacyl-tRNA biosynthesis, riboflavin metabolism, phenylpropanoid biosynthesis, and phenylpropanoid / flavonoid biosynthesis (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ea, b; Fig. S3; Tables S10, S11). These pathway enrichments suggest potential lineage-specific adaptations in secondary metabolite production in \u003cem\u003eS. stolonifera\u003c/em\u003e. Phylogenomic analysis using 799 single-copy orthologs confirmed that Saxifragaceae is sister to Grossulariaceae, with \u003cem\u003eS. stolonifera\u003c/em\u003e diverging from other Saxifragaceae species (\u003cem\u003eBergenia scopulosa\u003c/em\u003e, \u003cem\u003eTiarella polyphylla\u003c/em\u003e, and \u003cem\u003eChrysosplenium sinicum\u003c/em\u003e) approximately 59.39\u0026nbsp;million years ago (46.83\u0026ndash;71.99 MYA) (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ec). Gene family evolution analysis revealed 155 expanded and 17 contracted gene families in \u003cem\u003eS. stolonifera\u003c/em\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ec; Tables S12\u0026ndash;S14). Notably, expanded families were significantly enriched in flavonoid and phenylpropanoid biosynthesis pathways (Fig. S4a; Table S13). Anthocyanin accumulation not only reduces photodamage[\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e, \u003cspan citationid=\"CR62\" class=\"CitationRef\"\u003e62\u003c/span\u003e] but also underlies purple pigmentation in plants[\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e, \u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e]. Thus, the expansion of these gene families appears functionally linked to anthocyanin biosynthesis in \u003cem\u003eS. stolonifera\u003c/em\u003e.\u003c/p\u003e \u003cp\u003eWhole-genome duplications (WGDs) are recognized as key drivers of plant genome evolution[\u003cspan citationid=\"CR63\" class=\"CitationRef\"\u003e63\u003c/span\u003e, \u003cspan citationid=\"CR64\" class=\"CitationRef\"\u003e64\u003c/span\u003e]. Analysis of synonymous substitution rates (Ks) among paralogous gene pairs revealed two distinct peaks at ~\u0026thinsp;0.32 and ~\u0026thinsp;1.37, indicating two whole-genome duplication (WGD) events (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ed; Fig. S5a), suggesting two independent WGD events. The ancient peak (Ks\u0026thinsp;\u0026asymp;\u0026thinsp;1.37) corresponds to the core eudicot gamma whole-genome triplication, whereas the recent peak (Ks\u0026thinsp;\u0026asymp;\u0026thinsp;0.32) represents a lineage-specific WGD that occurred after the divergence of \u003cem\u003eS. stolonifera\u003c/em\u003e from other Saxifragaceae species (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ee; Figs. S5b, S6)).\u003c/p\u003e \u003cp\u003eMicrosynteny analyses further support two rounds of WGD along the \u003cem\u003eS. stolonifera\u003c/em\u003e lineage, with a 4:1 gene copy ratio observed between the \u003cem\u003eS. stolonifera\u003c/em\u003e and \u003cem\u003eVitis vinifera\u003c/em\u003e genomes (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ed, e). Higher chromosomal collinearity was observed between \u003cem\u003eS. stolonifera\u003c/em\u003e and the three Saxifragaceae species (\u003cem\u003eT. polyphylla\u003c/em\u003e, \u003cem\u003eB. scopulosa\u003c/em\u003e, and \u003cem\u003eC. sinicum\u003c/em\u003e), compared with \u003cem\u003eR. nigrum\u003c/em\u003e (Grossulariaceae) (Fig. S5b), reinforcing the earlier divergence of \u003cem\u003eS. stolonifera\u003c/em\u003e from \u003cem\u003eR. nigrum\u003c/em\u003e and its closer relationship with the other Saxifragaceae members. Comparative genomic analyses between \u003cem\u003eS. stolonifera\u003c/em\u003e and \u003cem\u003eB. scopulosa\u003c/em\u003e or \u003cem\u003eC. sinicum\u003c/em\u003e revealed syntenic depth ratios of 4:3 and 5:3, respectively (Fig. S6), indicating a shared WGD. Together, these findings provide strong evidence for the shared gamma-WGT event across Saxifragaceae and confirm an independent, lineage-specific WGD in \u003cem\u003eS. stolonifera\u003c/em\u003e following this ancient duplication.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cdiv id=\"Sec12\" class=\"Section2\"\u003e \u003ch2\u003eMorphological and physiological basis of leaf variegation\u003c/h2\u003e \u003cp\u003eTransverse sections of the green-leaf type (G), white-variegated type (W), and purple-variegated type (P) leaves were prepared (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e; Fig. S7\u0026ndash;S9). These sections revealed that air spaces were prevalent within the epidermal and palisade tissues along the central vein in both the white-variegated and purple-variegated types (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e, b, c, e, f, h, i, k, l; Fig. S7\u0026ndash;S8). Notably, in the white-variegated type, air spaces were observed among palisade cells (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eb, e, h, k; Fig. S7). In contrast, epidermal and palisade tissue cells in the green-leaf type were densely packed (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003ea, d, g, j; Fig. S9). The chlorophyll content did not differ significantly among the three leaf phenotypes, nor across different regions of individual leaves or developmental stages (Fig. S10; Table S15). Therefore, the appearance of white veins primarily results from air spaces located between the epidermal and palisade tissues along the central vein in the variegated types.\u003c/p\u003e \u003cp\u003eIn transverse sections, the spatial distribution of red pigments varied across leaf phenotypes. In the white-variegated type, pigment-containing cells were scattered throughout the mesophyll and hypoepidermal tissues (Fig. S7c), whereas in the green-leaf type, they were predominantly confined to the hypoepidermis (Fig. S9c). Neither the green-leaf nor the white-variegated type exhibited purple patches on the adaxial surface (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003ea\u0026ndash;c; Fig. S7a, S9a). Conversely, the purple-variegated type displayed accumulation of red pigments in the lower layers of the palisade tissue, leading to distinct purple patches visible on the upper leaf surface (Fig. S8a\u0026ndash;e). Collectively, the prominent white veins in the white-variegated type arise from intercellular air spaces between the epidermis and palisade tissue, coupled with altered cellular organization of the palisade layer. The conspicuous purple patches in the purple-variegated type result from the localized accumulation of red pigments in sub-palisade cells, while the associated white venation shares the same structural origin as in the white-variegated type. According to the classification of leaf variegations by Zhang \u003cem\u003eet al\u003c/em\u003e.[\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e], the white veins in the white-variegated type represent an air space-type variegation, whereas the purple-variegated type combines both pigment-based and air space-type mechanisms.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec13\" class=\"Section2\"\u003e \u003ch2\u003eMetabolomic profiling identifies anthocyanin glycosides as key pigments for purple-variegated type leaves\u003c/h2\u003e \u003cp\u003eTo identify the metabolites responsible for purple-variegated type, we profiled anthocyanins in leaves of the three phenotypes using LC-MS/MS. A total of 58 anthocyanin-related compounds were detected (Table S16). Comparative analysis revealed 10 differentially accumulated metabolites (DAMs) enriched in the anthocyanin and flavonol biosynthesis pathways (ko00942, ko00944), including cyanidin-3-O-glucoside, peonidin-3-O-galactoside, and quercetin-3-O-glucoside (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea; Table S17). Notably, these compounds accumulated at significantly higher levels in purple-variegated type compared to green-leaf type, white-variegated type leaves (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea; Table S17). All identified pigments are glycosylated forms of anthocyanidins, suggesting that glycosylation plays a critical role in their stabilization and accumulation. These results pinpoint anthocyanin glycosides as the chemical basis of the purple patches.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec14\" class=\"Section2\"\u003e \u003ch2\u003ePlant hormone profiling identifies for white-variegated type leaves\u003c/h2\u003e \u003cp\u003ePlant hormones profiling detected 46 hormones across eight classes (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eb; Table S18). Comparative analysis identified 35 differentially accumulated hormones, with the most pronounced differences observed between purple-variegated type and green-leaf type leaves (Table S19). Strikingly, both purple-variegated type and white-variegated type (W and P) accumulated significantly higher levels of multiple gibberellins (GA4, GA7, GA9, GA15, GA24) compared to green-leaf type leaves, in both young and first-stage adult leaves (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ec; Table S19). Among these, GA4 and GA7 are bioactive forms known to promote cell expansion[\u003cspan citationid=\"CR65\" class=\"CitationRef\"\u003e65\u003c/span\u003e]. The elevated GA levels in variegated leaves point to a potential role for gibberellin in mediating air space formation.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e \u003cb\u003eTranscriptomic analysis identifies\u003c/b\u003e \u003cb\u003eSsBZ1\u003c/b\u003e \u003cb\u003eand\u003c/b\u003e \u003cb\u003eSsGA3ox\u003c/b\u003e \u003cb\u003eas key regulators of pigment type and air space type leaf variegation\u003c/b\u003e\u003c/p\u003e \u003cp\u003eTo identify genes associated with the formation of air space-type and pigment-type variegation, RNA sequencing was performed on leaves across the three developmental stages of the three leaf phenotypes in \u003cem\u003eS. stolonifera\u003c/em\u003e. An average of 7.32 Gb clean reads per sample was obtained (Table S20), with an average mapping rate of 85%. A total of 18,923 differentially expressed genes (DEGs) were identified across all comparisons, with the highest number observed in the W2 vs. G2 comparison (6,209 DEGs; Fig. S11a; Table S21). Among these, 4,457 DEGs were commonly differentially expressed in at least two of the three stage specific comparisons between the purple-variegated and green-leaf types (YP-vs-YG, P1-vs-G1, and P2-vs-G2; Fig. S11b; Table S21). Similarly, 4,895 DEGs were shared across at least two comparisons between the white-variegated and green-leaf types (YW-vs-YG, W1-vs-G1, and W2-vs-G2; Fig. S11d; Table S21). Integrating these results with our high-quality \u003cem\u003eS. stolonifera\u003c/em\u003e genome, we focused on key regulatory genes involved in anthocyanin and gibberellin biosynthesis, reconstructing the anthocyanin (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e) and gibberellin (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003ea) pathways associated with leaf variegation.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eWe annotated 63 enzymatic genes involved in the anthocyanin biosynthesis pathway and visualized their expression patterns across three developmental stages of the three leaf phenotypes (Table S22). By reconstructing anthocyanin synthesis in \u003cem\u003eS. stolonifera\u003c/em\u003e leaves, we found that \u003cem\u003eSs4CL\u003c/em\u003e (\u003cem\u003eSst04G018150\u003c/em\u003e, \u003cem\u003eSst04G018170\u003c/em\u003e, \u003cem\u003eSst04G019880\u003c/em\u003e, \u003cem\u003eSst14G006060\u003c/em\u003e), \u003cem\u003eSsCHS\u003c/em\u003e (\u003cem\u003eSst12G013330, Sst08G008480\u003c/em\u003e), and \u003cem\u003eSsBZ1\u003c/em\u003e (\u003cem\u003eSst14G006050, Sst12G009310\u003c/em\u003e) were expressed at relatively high levels in fully developed leaves of the purple-variegated type (P2) (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e; Fig. S11 c; Table S22, S23). Notably, anthocyanidin 3-O-glucosyltransferase genes (\u003cem\u003eSsBZ1\u003c/em\u003e: \u003cem\u003eSst12G009310\u003c/em\u003e, \u003cem\u003eSst14G006050\u003c/em\u003e) showed elevated expression in both the white-variegated (W2) and purple-variegated (P2) leaf types, coinciding with the significant accumulation of key pigmentation metabolites: cyanidin-3,5-O-diglucoside, cyanidin-3-O-glucoside, peonidin-3-O-galactoside, peonidin-3-O-glucoside, peonidin-3,5-O-diglucoside, and quercetin-3-O-glucoside (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea; Table S17).\u003c/p\u003e \u003cp\u003eIn gibberellin biosynthesis, 24 enzymatically active genes were identified, and their expression profiles were analyzed across the three leaf phenotypes (Table S24). Reconstruction of the gibberellin pathway in \u003cem\u003eS. stolonifera\u003c/em\u003e revealed that gibberellin 3β-dioxygenase (\u003cem\u003eSsGA3ox\u003c/em\u003e: \u003cem\u003eSst05G017510\u003c/em\u003e) was highly expressed in all three phenotypes, particularly in the white-variegated type (groups W1 and W2), likely contributing to sustained high concentrations of bioactive GA4 (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003ea; Fig. S11 e; Table S25). These expression patterns are consistent with the significant accumulation of GA4, the primary active gibberellin, in these tissues (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eb, c, \u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003ea; Table S19). Collectively, these results suggest that the upregulated expression of GA 3-oxidases (\u003cem\u003eGA3ox\u003c/em\u003e), which maintains elevated levels of active gibberellin, is a key factor underlying white vein formation.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec15\" class=\"Section2\"\u003e \u003ch2\u003eWGCNA identifies a gibberellin-associated co-expression module with SsGA3ox\u003c/h2\u003e \u003cp\u003eTo gain further insights into the regulatory mechanisms underlying gibberellin content fluctuations during leaf variegation, weighted gene co-expression network analysis (WGCNA) was performed. A total of 9,960 DEGs were filtered and retained based on Fragments Per Kilobase of transcript per Million mapped reads (FPKM)\u0026thinsp;\u0026gt;\u0026thinsp;1 in more than 90% of the samples and showing significant differential expression across the three developmental stages of the three leaf phenotypes (Fig. S12). This analysis revealed 23 co-expression modules (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003eb, c, color-coded). Among these gene co-expression subnetworks (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003ec), several exhibited strong correlations with the synthesis and activation of the five gibberellins. Notably, the two bioactive gibberellins, GA4 and GA7, showed highly significant positive correlations with the salmon-coded module, with correlation coefficients of 0.9 and 0.94, respectively. Gene expression patterns within this module were characterized by high transcript levels in the white-variegated and purple-variegated leaf types (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003ea; Fig. S11 d, e). Within the salmon-coded module, the key gene \u003cem\u003eSsGA3ox\u003c/em\u003e (\u003cem\u003eSst05G017510\u003c/em\u003e) displayed significantly elevated expression in white-variegated leaves, particularly in fully expanded adult leaves (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003ea; Fig. S13). These findings support the hypothesis that upregulated \u003cem\u003eSsGA3ox\u003c/em\u003e expression maintains elevated levels of active gibberellins, thereby promoting rapid leaf expansion and contributing to the formation of intercellular air spaces between the epidermal and palisade tissues along the veins in air space type-leaf variegation.\u003c/p\u003e \u003c/div\u003e"},{"header":"Discussion","content":"\u003cp\u003eHere, we present a high-quality chromosome-scale genome assembly for \u003cem\u003eS. stolonifera\u003c/em\u003e. Ks distribution analyses provided compelling evidence for a shared gamma WGT event among Saxifragaceae species and revealed an independent WGD event in \u003cem\u003eS. stolonifera\u003c/em\u003e after the gamma-WGT. Comparative analysis of gene families across 25 species identified expanded gene families in \u003cem\u003eS. stolonifera\u003c/em\u003e that are highly enriched in phenylpropanoid biosynthesis, fatty acid metabolism, stilbenoid diarylheptanoid and gingerol biosynthesis, and flavonoid biosynthesis. Additionally, we found that \u003cem\u003eGA3ox\u003c/em\u003e genes in \u003cem\u003eS. stolonifera\u003c/em\u003e play a significant role in maintaining elevated levels of active gibberellins. Collectively, these findings suggest that high concentrations of active gibberellins may have been crucial for the development of air spaces between the epidermal and palisade tissues along leaf veins.\u003c/p\u003e \u003cp\u003e \u003cb\u003eLocalized anthocyanin accumulation determines the pigment type of\u003c/b\u003e \u003cb\u003eS. stolonifera\u003c/b\u003e\u003c/p\u003e \u003cp\u003eBioinformatic analysis of \u003cem\u003eS. stolonifera\u003c/em\u003e leaves indicated that the up-regulated expression of the anthocyanidin 3-O-glucosyltransferase gene (\u003cem\u003eSsBZ1\u003c/em\u003e: \u003cem\u003eSst12G009310\u003c/em\u003e, \u003cem\u003eSst14G006050\u003c/em\u003e) was a primary contributor to the formation of the purple-variegated phenotype. Anthocyanidin 3-O-glucosyltransferase is a key glycosyltransferase in plants[\u003cspan citationid=\"CR66\" class=\"CitationRef\"\u003e66\u003c/span\u003e]. Glycosylation plays a critical role in secondary metabolite biosynthesis by enhancing the stability, solubility, subcellular localization, and biological activity of conjugated compounds[\u003cspan additionalcitationids=\"CR67 CR68\" citationid=\"CR66\" class=\"CitationRef\"\u003e66\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR69\" class=\"CitationRef\"\u003e69\u003c/span\u003e]. Enhanced expression of this gene promotes the glucosylation of anthocyanidins into stable anthocyanins, thereby preserving their coloration[\u003cspan citationid=\"CR55\" class=\"CitationRef\"\u003e55\u003c/span\u003e]. Previous studies have shown that the anthocyanidin 3-O-glucosyltransferase gene is up-regulated in anthocyanin-rich, non-green leaves, such as those of \u003cem\u003eAlternanthera bettzickiana\u003c/em\u003e[\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e], \u003cem\u003eCamellia sinensis\u003c/em\u003e \u0026lsquo;Zijuan\u0026rsquo;[\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e], and \u003cem\u003ePerilla frutescens\u003c/em\u003e[\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e]. Similarly, in transcriptomic analyses of pigment-type variegation in \u003cem\u003eBegonia masoniana\u003c/em\u003e, the \u003cem\u003eUFGT\u003c/em\u003e gene (anthocyanidin 3-O-glucosyltransferase), involved in anthocyanin glycosylation is also up-regulated[\u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eIn this study, expression of \u003cem\u003eSsBZ1\u003c/em\u003e (\u003cem\u003eSst12G009310\u003c/em\u003e, \u003cem\u003eSst14G006050\u003c/em\u003e) was significantly up-regulated in the purple-variegated type compared with the green-leaf type (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e; Table S22). The content of visible anthocyanins was markedly higher in the purple-variegated type than in the green-leaf type, consistent with the elevated expression of \u003cem\u003eSsBZ1\u003c/em\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea, \u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003e; Fig. S11 c; Table S17). Therefore, up-regulation of anthocyanidin 3-O-glucosyltransferase likely underlies the formation of the purple-variegated phenotype.\u003c/p\u003e \u003cp\u003eFurthermore, compared with the white-variegated type, visible anthocyanins in the purple-variegated type were localized specifically within palisade tissue cells between the veins rather than scattered in the mesophyll cells (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e; Fig. S7\u0026ndash;9). This spatial pattern explains why both types contain similar levels of visible anthocyanins (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea), yet only the purple-variegated type exhibits distinct purple patches. The localized accumulation of anthocyanins in palisade cells constitutes the structural basis for the observed pigmentation pattern and has also been documented in other pigment-variegated plants, including \u003cem\u003eNephelaphyllum pulchrum\u003c/em\u003e, \u003cem\u003eStreptolirion volubile\u003c/em\u003e, and \u003cem\u003eVriesea saundersii\u003c/em\u003e[\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e]. In conclusion, up-regulation of the anthocyanidin 3-O-glucosyltransferase gene enhances anthocyanin accumulation, forming the molecular foundation of pigment-type variegation in \u003cem\u003eS. stolonifera\u003c/em\u003e. However, the molecular mechanisms governing the specific localization of anthocyanins in palisade tissue cells between veins remain unclear.\u003c/p\u003e \u003cdiv id=\"Sec17\" class=\"Section2\"\u003e \u003ch2\u003eHigh gibberellin levels promote air spaces formation\u003c/h2\u003e \u003cp\u003eBased on the chromosome-level genome sequence, combined with gene co-expression and plant hormone profiling, we reconstructed the gibberellin biosynthesis pathway and its regulatory network in \u003cem\u003eS. stolonifera\u003c/em\u003e. GAs are phytohormones that regulate multiple aspects of plant development[\u003cspan citationid=\"CR70\" class=\"CitationRef\"\u003e70\u003c/span\u003e], including leaf expansion[\u003cspan citationid=\"CR71\" class=\"CitationRef\"\u003e71\u003c/span\u003e], morphogenesis[\u003cspan citationid=\"CR72\" class=\"CitationRef\"\u003e72\u003c/span\u003e], and final leaf size[\u003cspan citationid=\"CR73\" class=\"CitationRef\"\u003e73\u003c/span\u003e]. Elevated GA concentrations accelerate cell division and leaf elongation rates[\u003cspan additionalcitationids=\"CR75\" citationid=\"CR74\" class=\"CitationRef\"\u003e74\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR76\" class=\"CitationRef\"\u003e76\u003c/span\u003e]. Hormone profiling revealed that increased levels of bioactive GAs, particularly GA4 and GA7, play a pivotal role in the formation of air space-type leaf variegation in \u003cem\u003eS. stolonifera\u003c/em\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ec). During early leaf development, significant accumulation of multiple gibberellins, including GA4, GA7, GA9, GA15, and GA24, was observed in both the purple-variegated and white-variegated types (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ec; Table S19).\u003c/p\u003e \u003cp\u003eGA3ox catalyzes the final step in gibberellin biosynthesis, converting precursor GAs into the bioactive forms GA1 and GA4[\u003cspan additionalcitationids=\"CR78 CR79 CR80 CR81\" citationid=\"CR77\" class=\"CitationRef\"\u003e77\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR82\" class=\"CitationRef\"\u003e82\u003c/span\u003e]. These enzymes are also known as gibberellin 3β-dioxygenases. Transcriptomic data demonstrated that the gibberellin 3β-dioxygenase gene (\u003cem\u003eSsGA3ox: Sst05G017510\u003c/em\u003e), which mediates this conversion, was highly expressed in white-variegated leaves, especially at the fully mature stage (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003ea; Fig. S11d e, S13). This sustained expression likely maintains elevated levels of active gibberellins, promoting rapid cell division and leaf expansion, key processes facilitating physical separation between epidermal and palisade tissue cells and resulting in air space formation along veins. Supporting this, WGCNA identified a co-expression module (salmon-colored) strongly correlated with GA4 and GA7 levels (Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e7\u003c/span\u003ec), and genes within this module, including \u003cem\u003eSsGA3ox\u003c/em\u003e, were highly expressed in variegated leaves (Fig. S12, S13). Together, these results indicate that gibberellin mediated cell expansion and tissue patterning form the structural basis of air space-type variegation. This mechanism enhances our understanding of leaf morphological diversity and provides a physiological and molecular explanation for the development of non-pigmented white veins in variegated leaves of \u003cem\u003eS. stolonifera\u003c/em\u003e.\u003c/p\u003e \u003cp\u003eIn sum, the high-quality reference genome of \u003cem\u003eS. stolonifera\u003c/em\u003e presented here serves as a valuable resource for studying evolutionary dynamics and phenotypic diversity within the genus \u003cem\u003eSaxifraga\u003c/em\u003e. Moreover, comprehensive transcriptomic, anthocyanin metabolic, and plant hormone profiling analyses provide novel insights into the molecular mechanisms underlying both pigment-type and air space-type variegation. The data generated in this study represent essential resources for future functional genomics and genetic investigations in \u003cem\u003eSaxifraga\u003c/em\u003e.\u003c/p\u003e \u003c/div\u003e"},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eEthics approval and consent to participate\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eNot applicable\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eConsent for publication\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eAll authors agree to the publication of this manuscript upon its acceptance.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eConflict of interest\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe authors declare no competing interests.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAvailability of data and materials\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe data that support the findings of this study have been deposited in the CNSA (accessed on 28 November 2024, https://db.cngb.org/cnsa/) of CNGBdb with accession code CNP0006535.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAcknowledgments\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eWe extend our sincere thanks to\u0026nbsp;Dr. Biao Xiong (Guizhou University) for guiding some data analyses.\u0026nbsp;Our gratitude goes to Quangang Xue, Jiayu Jin, and Boyang Xie (Shaoxing University) for their kind help with the materials planting. We thank TopEdit (www.topeditsci.com) for linguistic assistance during the preparation of this manuscript. We also appreciate the assistance of Dr. Bin Dong (Zhejiang A\u0026amp;F University) for his help in improving our manuscript.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAI declaration\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eDuring the preparation of this work, the author(s) used [DeepSeek] solely for language polishing and to improve readability. After using this tool, the author(s) reviewed and edited the content as needed and take(s) full responsibility for the content of the publication.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAuthor contributions\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eH.L., S.L., J.L., and J.Z. conceived and designed the project; H.L., J.Z., and S.L. provided financial support; S.L., J.L., and J.Z. collected the samples; S.L., J.L., and J.Z. performed the experiments; F.Z., H.X., F.T., J.Z., and S.L. analyzed the data; F.Z., H.X., F.T., and C.W. contributed to project discussion. H.L., J.Z., and S.L. wrote the manuscript. All authors contributed critically to the drafts and gave final approval for publication. All authors have read and agreed to the published version of the manuscript.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eFunding\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe project was supported by a grant from the National Natural Science Foundation of China (32300178) to S.L., the Startup Fund for Shaoxing University (13011001002/218) to J.Z., Biodiversity Investigation, Observation and Assessment Program (2019-2023) of Ministry of Ecology and Environment of China to H.L.\u003c/p\u003e\u003cp\u003e\u003cstrong\u003eSupplemental data\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eSupplementary data are available at \u003cem\u003eBMC Genomics\u0026nbsp;\u003c/em\u003eonline.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eZhang JH, Zeng JC, Wang XM, et al. A revised classification of leaf variegation types. Flora. 2020;272:151703.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang JH, Li JC, Zou L, et al. Transcriptome analysis of air space-type variegation formation in \u003cem\u003eTrifolium pratense\u003c/em\u003e. Int J Mol Sci. 2022;23:7794\u0026ndash;7794.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang JH. Classification of leaf variegation types in angiosperm and variegation related genes of \u003cem\u003eTrifolium pretense\u003c/em\u003e. East China Normal Univ. 2022.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLa Rocca N, Pupillo P, Puppi G et al. \u003cem\u003eErythronium dens-canis\u003c/em\u003e L. (Liliaceae): An unusual case of change of leaf mottling. \u003cem\u003ePlant Physiol Biochem.\u003c/em\u003e2014; 74: 108\u0026ndash;117.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSmith AP. Ecology of a leaf color polymorphism in a tropical forest species: habitat segregation and herbivory. Oecologia. 1986;69:283\u0026ndash;7.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSheue CR, Pao SH, Chien LF, et al. Natural foliar variegation without costs? The case of \u003cem\u003eBegonia\u003c/em\u003e. Ann Bot. 2012;109:1065\u0026ndash;74.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePao SH, Chesson P, Peng CI et al. Characteristics and ecological significance of natural foliar variegation from \u003cem\u003eBegonia formosana\u003c/em\u003e. In: The 99th Annual Meeting, 10th-15th (2014) August. Sacramento, California, America. Ecological Society of America.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShelef O, Summerfield L, Lev-Yadun S, et al. Thermal benefits from white variegation of \u003cem\u003eSilybum marianum\u003c/em\u003e leaves. Front Plant Sci. 2019;10:688.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi X, Kanakala S, He YH, et al. Physiological characterization and comparative transcriptome analysis of white and green leaves of \u003cem\u003eAnanas comosus\u003c/em\u003e var. \u003cem\u003ebracteatus\u003c/em\u003e. PLoS ONE. 2017;12:e0169838.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSun YH, Hung CY, Qiu J, et al. Accumulation of high OPDA level correlates with reduced ROS and elevated GSH benefiting white cell survival in variegated leaves. Sci Rep. 2017;7:1\u0026ndash;16.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang QM, Cui JG, Dai HY et al. Comparative transcriptome profiling of genes and pathways involved in leaf-patterning of \u003cem\u003eClivia miniata\u003c/em\u003e var. \u003cem\u003evariegata\u003c/em\u003e. \u003cem\u003eGene\u003c/em\u003e. 2018; 677: 280\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi WJ, Li HG, Shi LS, et al. Leaf color formation mechanisms in \u003cem\u003eAlternanthera bettzickiana\u003c/em\u003e elucidated by metabolite and transcriptome analyses. Planta. 2022;255:59.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eXie H, Zhu JY, Wang H, et al. An enhancer-transposable element from purple leaf tea varieties underlies the transition from evergreen to purple leaf color. Plant Commun. 2025;6(2):101176.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiu XN, Zhai YN, Liu JY, et al. Comparative transcriptome sequencing analysis to postulate the scheme of regulated leaf coloration in \u003cem\u003ePerilla frutescens\u003c/em\u003e. Plant Mol Biol. 2023;112:119\u0026ndash;42.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen J, Li Y, He D, et al. Cytological, physiological and transcriptomic analysis of variegated Leaves in \u003cem\u003ePrimulina pungentisepala\u003c/em\u003e offspring. BMC Plant Bio. 2022;22:419.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTanaka Y, Sasaki N, Ohmiya A. Biosynthesis of plant pigments: anthocyanins, betalains and carotenoids. Plant J. 2008;54:733\u0026ndash;49.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYang T, Xue ZZ, Li N. Reference genes selection and validation in \u003cem\u003eBegonia masoniana\u003c/em\u003e leaves of different developmental stages. Acta Hortic Sin. 2021;48:2251\u0026ndash;61.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhu J, Wang YZ, Wang QY, et al. The combination of DNA methylation and positive regulation of anthocyanin biosynthesis by MYB and bHLH transcription factors contributes to the petal blotch formation in Xibei tree peony. Hortic Res. 2023;10:uhad100.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAraguirang GE, Richter AS. Activation of anthocyanin biosynthesis in high light-what is the initial signal? New Phytol. 2022;236:2037\u0026ndash;43.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCui WH, Guan KY. Diversity of leaf variegation in Chinese Begonias. Plant Divers Resour. 2013;35:119\u0026ndash;27.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang H, Zhang P, Niu Y, et al. Genetic basis of camouflage in an alpine plant and its long-term co-evolution with an insect herbivore. Nat Ecol Evol. 2025;9:628\u0026ndash;38.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKong LY, Ou GJ, Wang CH. Application of \u003cem\u003eSaxifraga stolonifera\u003c/em\u003e Curt. to landscape architecture. Shaanxi For Sci Technol. 2014;2:77\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTang SM, Cai WQ, Zhang DM, et al. Effects of light intensites on morphological and physiological idexes of three ornamental cultivars of \u003cem\u003eSaxifraga stolonifera\u003c/em\u003e. Guihaia. 2023;43:699\u0026ndash;711.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi GF, He Y, Qin SH, et al. Optimization of chromosome preparation and karyotype analysis of polyploid watermelons. J South Agric. 2017;48:663\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMar\u0026ccedil;ais G, Kingsford C. A fast, lock-free approach for efficient parallel counting of occurrences of k-mers. Bioinformatics. 2011;27:764\u0026ndash;70.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiu B, Shi Y, Yuan J, et al. Estimation of genomic characteristics by analyzing k-mer frequency in de novo genome projects. arXiv. 2012. \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003e10.48550/arXiv.1308.2012.1308.2012v2\u003c/span\u003e\u003cspan address=\"10.48550/arXiv.1308.2012.1308.2012v2\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFlynn JM, Hubley R, Goubert C, et al. RepeatModeler2 for automated genomic discovery of transposable element families. PNAS. 2020;117:9451\u0026ndash;7.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBao Z, Eddy SR. Automated de novo identification of repeat sequence families in sequenced genomes. Genome Res. 2002;12:1269\u0026ndash;76.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePrice AL, Jones NC, Pevzner PA. De novo identification of repeat families in large genomes. Bioinformatics. 2005; Suppl 1: i351.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eEllinghaus D, Kurtz S, Willhoeft U. LTRharvest, an efficient and flexible software for de novo detection of LTR retrotransposons. BMC Bioinforma. 2008;9:18.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eXu Z, Wang H. LTR_FINDER: an efficient tool for the prediction of full-length LTR retrotransposons. Nucleic Acids Res. 2007;35:W265\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTarailo-Graovac M, Chen NS. Using RepeatMasker to identify repetitive elements in genomic sequences. Curr Protocols Bioinf 2009; 4: 4.10. 11\u0026ndash;14.10. 14.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBenson G. Tandem repeats finder: a program to analyze DNA sequences. Nucleic Acids Res. 1999;27:573\u0026ndash;80.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBeier S, Thiel T, M\u0026uuml;nch T, et al. Misa-web: a web server for microsatellite prediction. Bioinformatics. 2017;33:2583\u0026ndash;5.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCantarel BL, Korf I, Robb SMC. MAKER: An easy-to-use annotation pipeline designed for emerging model organism genomes. Genome Res. 2008;18:188\u0026ndash;96.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKeilwagen J, Wenk M, Erickson JL, et al. Using intron position conservation for homology-based gene prediction. Nucleic acids res. 2016;44:e89\u0026ndash;89.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eStanke M, Keller O, Gunduz I, et al. AUGUSTUS: ab initio prediction of alternative transcripts. Nucleic Acids Res. 2006;34:W435\u0026ndash;9.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKorf I. Gene finding in novel genomes. BMC Bioinforma. 2004;5:59.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKim D, Langmead B, Salzberg SL. HISAT: A fast spliced aligner with low memory requirements. Nat Methods. 2015;12:357\u0026ndash;60.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePertea M, Kim D, Pertea GM, et al. Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat Protoc. 2016;11:1650\u0026ndash;67.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTang S, Lomsadze A, Borodovsky M. Identification of protein coding regions in RNA transcripts. Nucleic Acids Res. 2015;43:e78.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHaas BJ, Delcher AL, Mount SM, et al. Improving the Arabidopsis genome annotation using maximal transcript alignment assemblies. Nucleic Acids Res. 2003;31:5654\u0026ndash;66.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNip KM, Chiu R, Yang C, et al. RNA-Bloom enables reference-free and reference-guided sequence assembly for single-cell transcriptomes. Genome Res. 2020;30:1191\u0026ndash;200.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHaas BJ, Salzberg SL, Zhu W, et al. Automated eukaryotic gene structure annotation using EVidenceModeler and the Program to Assemble Spliced Alignments. Genome Biol. 2008;9:R7.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSim\u0026atilde;o FA, Waterhouse RM, Ioannidis P, et al. BUSCO: assessing genome assembly and annotation completeness with single-copy orthologs. Bioinformatics. 2015;31:3210\u0026ndash;2.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eEmms DM, Kelly S. OrthoFinder: phylogenetic orthology inference for comparative genomics. Genome boil. 2019;20:1\u0026ndash;14.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMistry J, Chuguransky S, Williams L, et al. Pfam: The protein families database in 2021. Nucleic Acids Res. 2021;49:D412\u0026ndash;9.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKatoh K, Asimenos G, Toh H. Multiple alignment of DNA sequences with MAFFT. Methods Mol Biol. 2009;537:39\u0026ndash;64.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNguyen LT, Schmidt HA, von Haeseler A, et al. IQ-TREE: a fast and effective stochastic algorithm for estimating maximum-likelihood phylogenies. Mol Biol Evol. 2015;32:268\u0026ndash;74.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHan MV, Thomas GWC, Lugo-Martinez J, et al. Estimating gene gain and loss rates in the presence of error in genome assembly and annotation using CAFE 3. Mol Biol Evol. 2013;30:1987\u0026ndash;97.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSun PC, Jiao BB, Yang YZ, et al. A user-friendly toolkit for evolutionary analyses of whole-genome duplications and ancestral karyotypes. Mol Plant. 2022;15:1841\u0026ndash;51.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMeng Y, Zhang H, Fan Y, et al. Anthocyanins accumulation analysis of correlated genes by metabolome and transcriptome in green and purple peppers (\u003cem\u003eCapsicum annuum\u003c/em\u003e). BMC Plant Biol. 2022;22:358.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang Q, Zhu J, Li B, et al. Functional identification of anthocyanin glucosyltransferase genes: a Ps3GT catalyzes pelargonidin to pelargonidin 3-O-glucoside painting the vivid red flower color of Paeonia. Planta. 2023;257:56.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGeng Z, Liu M, Wang Y, et al. Transcriptomic and metabolomic analyses reveal CmMYB308 as a key regulator in the pink flower color variation of \u0026lsquo;Dante Purple\u0026rsquo; chrysanthemum. Plant Cell Rep. 2024;43:157.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eXia XM, Du HL, Hu XD, et al. Genomic insights into adaptive evolution of the species-rich cosmopolitan plant genus \u003cem\u003eRhododendron\u003c/em\u003e. Cell Rep. 2024;43:114745.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiu FH, Zhao Y, Wang XR, et al. Regulatory roles of GA and sugar metabolism genes with metabolomic profiling during seed dormancy release in \u003cem\u003eGleditsia sinensis\u003c/em\u003e. BMC Plant Biol. 2025;25:1231.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen C, Wu Y, Li J, et al. TBtools-II: a one for all, all for one bioinformatics platform for biological big-data mining. Mol Plant. 2023;16:1733\u0026ndash;42.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLove MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang X, Liang S, Yang WQ, et al. MetMiner: A user-friendly pipeline for large-scale plant metabolomics data analysis. J Integr Plant Biol. 2024;66:2329\u0026ndash;45.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBendich AJ, Rogers SO. The biological and evolutionary consequences of competition between DNA sequences that benefit the cell and DNA sequences that benefit themselves. Nucleic Acids Res. 2025;53:gkaf589.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang F, Jiang Z, Gao J, et al. Analysis of the \u003cem\u003eRehmannia chingii\u003c/em\u003e geneome identifies RcCYP72H7 as an epoxidase in iridoid glycoside biosynthesis. Nat Commun. 2025;16:6035.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSolovchenko AE, Merzlyak MN. Screening of visible and UV radiation as a photoprotective mechanism in plants. Russ J Plant Physiol. 2008;55:719\u0026ndash;37.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eJiao YN, Wickett NJ, Ayyampalayam S, et al. Ancestral polyploidy in seed plants and angiosperms. Nature. 2011;473:97\u0026ndash;100.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAdams K. Genomic clues to the ancestral flowering plant. Science. 2013;342:1456\u0026ndash;7.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHauvermale AL, Ariizumi T, Steber CM. Gibberellin signaling: a theme and variations on DELLA repression. Plant physiol. 2012;160:83\u0026ndash;92.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYin RH, Messner B, Faus-Kessler T, et al. Feedback inhibition of the general phenylpropanoid and flavonol biosynthetic pathways upon a compromised flavonol-3-O-glycosylation. J Exp Bot. 2012;63:2465\u0026ndash;78.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi Y, Baldauf S, Lim EK, et al. Phylogenetic analysis of the UDP-glycosyltransferase multigene family of \u003cem\u003eArabidopsis thaliana\u003c/em\u003e. J Biol Chem. 2001;276:4338\u0026ndash;43.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMe\u0026szlig;ner B, Thulke O, Sch\u0026auml;ffner AR. Arabidopsis glucosyltransferases with activities toward both endogenous and xenobiotic substrates. Planta. 2003;217:138\u0026ndash;46.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBowles D, Lim EK, Poppenberger B, et al. Glycosyltransferases of lipophilic small molecules. Annu Rev Plant Biol. 2006;57:567\u0026ndash;97.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSalazar-Cerezo S, Mart\u0026iacute;nez-Montiel N, Garc\u0026iacute;a-S\u0026aacute;nchez J, et al. Gibberellin biosynthesis and metabolism: A convergent route for plants, fungi and bacteria. Microbiol Res. 2018;208:85\u0026ndash;98.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi J, Gao HW, Jiang JS, et al. Overexpression of the \u003cem\u003eGalega orientalis\u003c/em\u003e gibberellin receptor improves biomass production in transgenic tobacco. Plant Physiol Biochem. 2013;73:1\u0026ndash;6.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGao LW, Lyu SW, Tang J, et al. Genome-wide analysis of auxin transport genes identifies the hormone responsive patterns associated with leafy head formation in Chinese cabbage. Sci Rep. 2017;7:42229.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMu XH, Chen QW, Wu XY, et al. Gibberellins synthesis is involved in the reduction of cell flux and elemental growth rate in maize leaf under low nitrogen supply. Environ Exp Bot. 2018;150:198\u0026ndash;208.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWenzel CL, Williamson RE, Wasteneys GO, et al. Gibberellin-induced changes in growth anisotropy precede gibberellin-dependent changes in cortical microtubule orientation in developing epidermal cells of barley leaves. Kinematic and cytological studies on a gibberellin-responsive dwarf mutant, M489. Plant Physiol. 2000;124:813\u0026ndash;22.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNelissen H, Rymen B, Jikumaru Y, et al. A local maximum in gibberellin levels regulates maize leaf growth by spatial control of cell division. Curr Biol. 2012;22:1183\u0026ndash;7.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang XJ, He LL, Zhao BL, et al. \u003cem\u003eDwarf and Increased Branching 1\u003c/em\u003e controls plant height and axillary bud outgrowth in \u003cem\u003eMedicago truncatula\u003c/em\u003e. J Exper Bot. 2020;71:6355\u0026ndash;65.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLester DR, Ross JJ, Davies PJ, et al. Mendel\u0026rsquo;s stem length gene (Le) encodes a gibberellin 3 beta-hydroxylase. Plant Cell. 1997;9:1435\u0026ndash;43.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eItoh H, Ueguchi-Tanaka M, Sentoku N, et al. Cloning and functional analysis of two gibberellin 3β-hydroxylase genes that are differently expressed during the growth of rice. PNAS. 2001;98:8909\u0026ndash;14.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDalmadi A, Kal\u0026oacute; P, Jakab J, et al. Dwarf plants of diploid \u003cem\u003eMedicago sativa\u003c/em\u003e carry a mutation in the gibberellin 3-beta-hydroxylase gene. Plant Cell Rep. 2008;27:1271\u0026ndash;9.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHedden P, Thomas SG. Gibberellin biosynthesis and its regulation. Biochem J. 2012;444:11\u0026ndash;25.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen Y, Hou M, Liu L, et al. The maize DWARF1 encodes a gibberellin 3-oxidase and is dual localized to the nucleus and cytosol. Plant Physiol. 2014;166:2028\u0026ndash;39.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePhillips AL, Huttly AK, Alarc\u0026oacute;n-Reverte R, et al. \u003cem\u003eGIBBERELLIN 3-OXIDASE\u003c/em\u003e genes regulate height and grain size in bread wheat. J Exp Bot. 2025;76:3345\u0026ndash;58.\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":false,"highlight":"","institution":"","isAcceptedByJournal":false,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"bmc-genomics","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":false,"externalIdentity":"gics","sideBox":"Learn more about [BMC Genomics](http://bmcgenomics.biomedcentral.com/)","snPcode":"","submissionUrl":"https://www.editorialmanager.com/gics","title":"BMC Genomics","twitterHandle":"#BMCGenomics","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"em","reportingPortfolio":"BMC Series","inReviewEnabled":true,"inReviewRevisionsEnabled":true},"keywords":"leaf variegation, anatomy, Saxifraga stolonifera, gibberellins, transcriptomics","lastPublishedDoi":"10.21203/rs.3.rs-9460163/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-9460163/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003e \u003cem\u003eSaxifraga stolonifera\u003c/em\u003e Curtis is a shade-tolerant, variegated-leaf herb with ornamental and medicinal value that naturally displays three stable leaf phenotypes: green, white-variegated, and purple-variegated. These phenotypes simultaneously exhibit both air space-type and pigment-type variegations. However, the key pathways and regulatory networks underlying their formation remain largely unknown. Here, we generated a 2.01 Gb high-quality chromosome-scale genome assembly (2n\u0026thinsp;=\u0026thinsp;36). Comparative genomic analyses revealed a recent whole-genome duplication specific to \u003cem\u003eS. stolonifera\u003c/em\u003e following the shared γ-triplication event, along with 155 expanded gene families enriched in flavonoid and phenylpropanoid biosynthesis pathways. Metabolomic profiling identified 58 anthocyanin-related compounds, among which 10 key pigments including cyanidin-3-O-glucoside, peonidin-3-O-galactoside, and quercetin-3-O-glucoside, were responsible for the purple patches. Their accumulation corresponded with the up-regulation of the anthocyanidin 3-O-glucosyltransferase gene \u003cem\u003eSsBZ1\u003c/em\u003e (\u003cem\u003eSst12G009310\u003c/em\u003e, \u003cem\u003eSst14G006050\u003c/em\u003e). Hormone and transcriptome analyses showed that white-variegated leaves accumulate high levels of bioactive gibberellins GA4 and GA7, driven by increased expression of gibberellin 3β-dioxygenase \u003cem\u003eSsGA3ox\u003c/em\u003e (\u003cem\u003eSst05G017510\u003c/em\u003e). Weighted gene co-expression network analysis (WGCNA) further identified a GA-associated module (salmon) with \u003cem\u003eSsGA3ox\u003c/em\u003e as the hub gene, promoting cell expansion and generating air spaces between epidermal and palisade tissues along the veins. Collectively, our high-quality genome, metabolome, and transcriptome resources demonstrate that \u003cem\u003eSsGA3ox\u003c/em\u003e-mediated gibberellin biosynthesis drives air space-type leaf variegation, whereas \u003cem\u003eSsBZ1\u003c/em\u003e-controlled anthocyanin glycosylation produces pigment-type leaf variegation. These findings provide an integrative omics framework for dissecting leaf variegation mechanisms in \u003cem\u003eSaxifraga\u003c/em\u003e and other ornamental plants.\u003c/p\u003e","manuscriptTitle":"Multi-omics reveals divergent regulation of anthocyanin glycosylation and gibberellin biosynthesis underlying leaf variegation in Saxifraga stolonifera","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2026-05-11 19:40:24","doi":"10.21203/rs.3.rs-9460163/v1","editorialEvents":[{"type":"communityComments","content":0},{"type":"reviewerAgreed","content":"273289778479265714058125319122522357430","date":"2026-05-16T02:56:56+00:00","index":"hide","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2026-05-10T18:07:29+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"57838931435626891902881849608034108448","date":"2026-05-04T16:33:54+00:00","index":"hide","fulltext":""},{"type":"reviewersInvited","content":"","date":"2026-05-04T02:51:59+00:00","index":"","fulltext":""},{"type":"editorAssigned","content":"","date":"2026-05-04T02:50:12+00:00","index":"","fulltext":""},{"type":"editorInvited","content":"","date":"2026-04-29T19:03:37+00:00","index":"","fulltext":""},{"type":"checksComplete","content":"","date":"2026-04-29T16:22:00+00:00","index":"","fulltext":""},{"type":"submitted","content":"BMC Genomics","date":"2026-04-29T14:45:36+00:00","index":"","fulltext":""}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"bmc-genomics","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":false,"externalIdentity":"gics","sideBox":"Learn more about [BMC Genomics](http://bmcgenomics.biomedcentral.com/)","snPcode":"","submissionUrl":"https://www.editorialmanager.com/gics","title":"BMC Genomics","twitterHandle":"#BMCGenomics","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"em","reportingPortfolio":"BMC Series","inReviewEnabled":true,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"821beac0-ada4-490c-81a2-5baff583af64","owner":[],"postedDate":"May 11th, 2026","published":true,"recentEditorialEvents":[{"type":"reviewerAgreed","content":"273289778479265714058125319122522357430","date":"2026-05-16T02:56:56+00:00","index":70,"fulltext":""},{"type":"editorInvitedReview","content":"","date":"2026-05-10T18:07:29+00:00","index":49,"fulltext":""},{"type":"reviewerAgreed","content":"57838931435626891902881849608034108448","date":"2026-05-04T16:33:54+00:00","index":45,"fulltext":""},{"type":"reviewersInvited","content":"40","date":"2026-05-04T02:51:59+00:00","index":"","fulltext":""},{"type":"editorAssigned","content":"","date":"2026-05-04T02:50:12+00:00","index":"","fulltext":""},{"type":"editorInvited","content":"","date":"2026-04-29T19:03:37+00:00","index":"","fulltext":""},{"type":"checksComplete","content":"","date":"2026-04-29T16:22:00+00:00","index":"","fulltext":""},{"type":"submitted","content":"BMC Genomics","date":"2026-04-29T14:45:36+00:00","index":"","fulltext":""}],"rejectedJournal":[],"revision":"","amendment":"","status":"under-review","subjectAreas":[],"tags":[],"updatedAt":"2026-05-11T19:40:24+00:00","versionOfRecord":[],"versionCreatedAt":"2026-05-11 19:40:24","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-9460163","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-9460163","identity":"rs-9460163","version":["v1"]},"buildId":"XKTyCvWXoU3ODBz1xrDgd","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}