# Fade into you: genetic control of pigmentation patterns in red-flesh apple (Malus domestica)

Pierre Bouillon, Etienne Belin, Anne-Laure Fanciullino, Sandrine Balzergue, Sylvain Hanteville, Yao Letekoma, Maryline Cournol, Fatima Faris, Andréa Bouanich, Dimitri Bréard, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

Pierre Bouillon, Etienne Belin, Anne-Laure Fanciullino, Sandrine Balzergue, Sylvain Hanteville, et al. Fade into you: genetic control of pigmentation patterns in red-flesh apple (Malus domestica). Frontiers in Plant Science, 2025, 15, pp.1462545. 10.3389/fpls.2024.1462545. hal-05081217

HAL Id: hal-05081217

https://hal.inrae.fr/hal-05081217v1

Submitted on 23 May 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

Check for updates

OPEN ACCESS

EDITED BY

Boas Pucker,
Technical University of Braunschweig,
Germany

REVIEWED BY

Saloni Sharma,
University of Granada, Spain
Asif Bashir Shikari,
Sher-e-Kashmir University of Agricultural
Sciences and Technology of Kashmir, India
Zhanwu Dai,
Chinese Academy of Sciences (CAS), China
Laura Dougherty,
Cornell University, United States

*CORRESPONDENCE

Jean-Marc Celton

<figure><img src="image_4.png" /></figure>

jean-marc.celton@inrae.fr

RECEIVED 10 July 2024

ACCEPTED 03 December 2024

PUBLISHED 13 January 2025

CITATION

Bouillon P, Belin E, Fanciullino A-L, Balzergue S, Hanteville S, Letekoma Y, Cournol M, Faris F, Bouanich A, Bréard D, Bernard F and Celton J-M (2025) Fade into you: genetic control of pigmentation patterns in red-flesh apple (Malus domestica). *Front. Plant Sci.* 15:1462545.
doi: 10.3389/fpls.2024.1462545

COPYRIGHT

© 2025 Bouillon, Belin, Fanciullino, Balzergue, Hanteville, Letekoma, Cournol, Faris, Bouanich, Bréard, Bernard and Celton. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Fade into you: genetic control of pigmentation patterns in red-flesh apple (Malus domestica)

Pierre Bouillon¹,², Etienne Belin¹, Anne-Laure Fanciullino¹, Sandrine Balzergue¹,³, Sylvain Hanteville¹, Yao Letekoma¹, Maryline Cournol¹, Fatima Faris¹, Andréa Bouanich¹, Dimitri Bréard⁴, Frédéric Bernard² and Jean-Marc Celton¹*

¹Univ Angers, Institut Agro, INRAE, IRHS, SFR QUASAV, Angers, France, ²IFO, Seiches sur le Loir, France, ³Analyses des Acides Nucléiques (ANAN), SFR QUASAV, Angers, France, ⁴Univ Angers, Substances d'Origine Naturelle et Analogues Structuraux (SONAS), SFR QUASAV, Angers, France

The genetic basis of type 1 red-flesh color development in apple (Malus domestica) depends upon a particular allele of the *MdMYB10* gene. Interestingly, type 1 red-flesh apples are fully red after fruit set, but anthocyanin pigmentation in apple fruit cortex may decrease during fruit growth and maturation, leading to variable red patterning and intensities in the mature cortical flesh. We developed a histogram-based color analysis method to quantitatively estimate pigmentation patterns. This methodology was applied to investigate the phenotypic diversity in four hybrid F1 families segregating for red-flesh color. Pigmentation patterns were found to be heritable allowing the identification of a new locus by QTL analysis. To further investigate the mechanisms involved in the spatial deposition of anthocyanin, metabolome, transcriptome and methylome comparisons between white and red flesh areas within the red-flesh genotype cv. 'R201' exhibiting flesh pigmentation patterns, was performed. Wide-targeted analysis showed that white-flesh areas accumulate more dihydrochalcones and hydroxycinnamic acids than red-flesh areas while red-flesh areas accumulate more flavonoids. Anthocyanin biosynthesis genes and anthocyanin positive regulators (MBW complex) were up-regulated in red-flesh areas, while a reduction in anthocyanin storage, transport and stability (increase of pH, down-regulation of *MdGSTU22*) and an increase in phenolic catabolism were concomitant with color fading process in white-flesh areas. Expression of *MdGSTU22* was linked to a differentially methylated region (DMR) suggesting a potential environmental effect on the epigenetic control of gene expression involved in color fading. Altogether, these results provide the first characterization and functional identification of color fading in apple fruit flesh.

KEYWORDS

red-flesh apple, *Malus domestica*, anthocyanin, pigmentation pattern, apple genetics
# 1 Introduction

Plant pigments have long fascinated scientists who have studied their chemical and physical properties, their biosynthesis, and their physiological and ecological functions (Grotewold, 2006). Anthocyanins, a class of flavonoids, are water-soluble phenolic compounds and pigments that provide vibrant hues in the plant kingdom, from orange/red to violet/blue (Tanaka et al., 2008). The sequential steps leading to biosynthesis of anthocyanins have been well described and are highly conserved among species (Buhrman et al., 2022). Anthocyanins synthesis begins with the condensation of three molecules of malonyl-CoA and one of p-coumaroyl-CoA from phenylalanine via the general phenylpropanoid pathway. Malonyl-CoA and p-coumaroyl-CoA are condensed by a chalcone synthase (CHS). Four enzymes in flavonoid biosynthesis, chalcone isomerase (CHI), flavanone 3-hydroxylase (F3H), dihydroflavonol 4-reductase (DFR) and anthocyanidin synthase (ANS or LDOX), are then required to achieve anthocyanidin synthesis (Zhang et al., 2014). Stable anthocyanins are then glycosylated on the 3-position by a (UDP)-glucose: flavonoid 3-O-glucosyltransferase (UFGT) (Liu et al., 2021).

Following primary synthesis, anthocyanins are transported to the vacuole where they accumulate and participate in the final color acquisition. Glutathione S-transferases (GSTs) (Zhao et al., 2021b), ATP-binding cassette (ABC-transporters) (Xiang et al., 2024) and toxic compound extrusion MATE-transporters (Zhang et al., 2021) are key enzymes that enable anthocyanin translocation from cytoplasm into the vacuole (Buhrman et al., 2022). In peach fruit (Prunus persica), loss-of-function alleles of the glutathione S-transferase (GST) gene causes anthocyanin deficiency in flower and fruit skin (Lu et al., 2021). In apple, MdGSTU12 was shown to play an important role in the regulation of anthocyanin accumulation in fruit skin (Zhao et al., 2021b). In kiwifruit, methylation of AcGST1 is associated with differential accumulation between the inner and outer pericarp (Zhang et al., 2023). Indeed, a decrease in the methylation level of the AcGST1 promoter may contribute to accumulation of anthocyanin in the outer pericarp of the fruit.

In addition to structural genes involved in anthocyanin biosynthetic pathway, several transcription factors (TFs) have been found to regulate anthocyanin accumulation (Allan et al., 2008). The ternary MBW protein complex comprising MYB, bHLH, and WD40 proteins is considered to be a major transcriptional regulator of anthocyanin synthesis (Chen et al., 2021). R2R3-MYB activators can bind, directly or within the MBW complex, to the promoters of anthocyanin biosynthetic genes (ABGs), thereby increasing their expression and modulating subsequent anthocyanin accumulation (Yan et al., 2021). In contrast, a large number of anthocyanin repressors have been identified in the last decade (LaFountain and Yuan, 2021). Some R2R3-MYB and R3 repressors are generally involved in interfering with MBW complex to inhibit activation of the anthocyanin regulatory network (LaFountain and Yuan, 2021), and activator-repressor MYB pairs have been found to control the formation of heterogeneous pigmentation region in many angiosperm species (Davies et al., 2012).

The color fading process, described as the loss of anthocyanin pigmentation in plant tissues and organs at particular

developmental stages, is well-known in ornamental species, particularly in flowers (Han et al., 2020; Huang et al., 2024). While regulatory mechanisms of anthocyanin biosynthesis are well described in many plant species, the genetic and biochemical mechanisms that lead to anthocyanin degradation are less well understood. In *Malus halliana*, repression of *MhMYB10* was associated with flower color fading (Han et al., 2020), suggesting that the increased methylation level of MYB10 slowed the anthocyanin accumulation process by repressing its expression. The importance of micro-RNAs (miRNAs) in anthocyanin regulation has also been characterized in *Malus Crabapple* where numerous miRNAs have been linked to petal color fading (Rong et al., 2023). In pear leaves and flowers, overexpression of a laccase through MYB TF activation resulted in color fading (Zhao et al., 2021a).

Study of pigmentation patterning is an exciting area (Fairnie et al., 2022) that can be traced back to the Reaction-Diffusion model established by Alan Turing (Kondo and Miura, 2010) that finds echoes in recent genetic models proposed on Monkeyflower (Mimulus) (Yuan, 2018; Ding et al., 2020). Color patterning can be dependant upon an activator-repressor system involving in the upregulation of an inhibitor during later steps of tissue development. For example, an activator-repressor system was observed in red-flesh kiwifruit (*Actinidia* spp.) where an interplay between positive regulation of MYB TF and negative regulation by mircoRNAs (miRNAs) regulated anthocyanin synthesis and distribution (Wang et al., 2022). An activator-repressor system was also identified in peach where a R2R3-MYB repressor, PpMYB18, competes with PpMYB10 for binding to bHLH in a fine-tuning regulatory loop to balance anthocyanin accumulation (Zhou et al., 2018). In contrast, traits without repeating pattern are generally explained in terms of positional specification, in which overlapping zones of gene expression regulate pigmentation pattern (Zheng et al., 2021). While biological functions of flower pigmentation pattern are well-studied in pollinators attraction (Moyers et al., 2017; Todesco et al., 2022) and defense against herbivore (Lev-Yadun, 2016), the biological process that leads to pigmentation pattern in fruits are much less understood. Recent progress in image analysis opens up a new field of study for pigmentation pattern with recent applications in fruit species such as strawberry (*Fragaria × ananassa*) (Zingaretti et al., 2021; Denoyes et al., 2023). In apple, fruit color is associated with anthocyanin accumulation (Allan et al., 2008). The major apple phenolic compounds are flavonols, monomeric and oligomeric flavan-3-ols, dihydrochalcones, hydroxycinnamic acids and anthocyanins, with cyanidin-3-galactoside accounting for about 80% of the total amount of anthocyanins in apple (Chen et al., 2021). Some genotypes display red pigmentation in the fleshy part of the fruit leading to the 'red-flesh' trait. Two different types of red-flesh apples have been characterized. The type 1 red-flesh phenotype accumulates anthocyanins throughout the fruit, from fruit set through maturity and is associated with red pigmentation of leaves, stems, roots, and flowers. The type 2 red-flesh apple phenotype is characterized by green vegetative tissues, yellow-orange fruit skin, and a red pigmentation that occurs only in the fruit flesh during later stages of fruit development (Chagné et al., 2012). Type 1 red-flesh is dependent upon the presence of a particular allele of the *MdMYB10*
gene (historically termed “R” locus) (Chagné et al., 2007; Espley et al., 2007). This gene, characterized by a minisatellite-like structure, is required for red-flesh pigmentation by upregulating anthocyanin late biosynthetic genes (Espley et al., 2009). Other genetic factors have been identified through recent genetic analyses (Kumar et al., 2022). Transcriptomic analyses identified differentially expressed genes (DEGs) associated with the flavonoid pathway between red-flesh and white-flesh apples (Wang et al., 2015). Among them, the WRKY-family TF MdWRKY11 seems to be a major regulator of flavonoid and anthocyanin synthesis (Wang et al., 2018), and is involved in vacuolar acidity (Verweij et al., 2016; Lloyd et al., 2017). Recently, Li et al (Li et al., 2022), observed that ABGs were more abundantly expressed in the red areas of the fruit than in the white sections. They identified, through transcriptomic analysis, DEGs related to phytohormones and involved in anthocyanin biosynthesis regulation, and proposed that some auxin and ethylene response TFs could participate in anthocyanin biosynthesis. MdNAC42 was observed to interact with MdMYB10 in red-flesh apple, and up-regulated the expression of several apple ABGs, leading to an increased fruit coloration (Zhang et al., 2020). MdMYB9 also participates in anthocyanin and proanthocyanidin synthesis through a HY5-miR858-MYB loop (Li et al., 2023). MYB repressors of anthocyanin synthesis have been successfully identified in apple. MrMYB44 caused color fading of Malus radiant leaves (Meng et al., 2022), and was closely linked to fruit malate accumulation (Jia et al., 2021). Interestingly, in strawberry, FaMYB44 was found to interact with FaMYB10 during development and regulate sugar/acid accumulation (Wei et al., 2018). In type 1 red-flesh apple, MdMYB6, MdMYB16, MdMYB17 and MdMYB111 are responsible for downregulation of anthocyanin structural genes (Xu et al., 2020; Wang and Chen, 2021).

Red-flesh fruits are fully red after fruit set (Espley et al., 2012), however, it has been observed that anthocyanin pigmentation in apple fruit cortex may decrease during the course of fruit growth and maturation, leading to variable red patterning and intensities in the cortical flesh (van Nocker et al., 2011; Volz et al., 2014; Bouillon et al., 2024a). This study aimed to identify genomic regions associated with flesh pigmentation pattern in apple. Red-flesh pigmentation pattern was found to be heritable in a hybrid progeny, underlying a potential genetic control of this trait (Bouillon et al., 2024a). In this study, a QTL analysis approach was combined with a differential transcriptome and a methylome analysis to decipher the genetic basis of pigmentation pattern in red-flesh apple fruit. A pedigree-based QTL analysis (PBA-QTL) was conducted using four F1 red-flesh apple progenies (256 individuals), supplemented with genetic data from parents, grandparents and founders (totaling 347 individuals) to characterize the genetic architecture of flesh pigmentation pattern. To conduct this analysis, an innovative phenotyping method was developed and implemented to quantitatively estimate flesh pigmentation pattern in red-flesh fruits. Comparison between inner white-flesh (WF) and outer red-flesh (RF) metabolome, transcriptome and methylome of a red-flesh cultivar (cv. 'R201') that exhibited flesh pigmentation pattern was then performed.

# 2 Material and methods

## 2.1 Plant material

A pedigree-based QTL analysis (PBA-QTL) was performed on 256 individuals from four interconnected (common parent and/or grandparent) F1 families (Bouillon et al., 2024a), supplemented with genetic data from parents, grandparents and founders (totaling 347 individuals). Three families segregated for red-flesh (RF1-1 to RF1-3, plantation year: 2018) and one family was a cross between two red-flesh parents (RF1-5, plantation year: 2017). Hybrids were selected at the seedling stage for red-leaves color [phenotypic marker of MdMYB10 (Chagné et al., 2007)] and, consequently, possessed at least one R6-MdMYB10 allele. Fruit harvest was conducted for three years (2021, 2022 and 2023) from August to end-October in IFO orchard (L'Anguicherie, 49140 Seiches-sur-Le Loir, France/GCS: 47°37'52.5"N 0°19'38.4"W). For each genotype harvested at maturity (brix values varying from 13 to 22; starch index between 6 and 8), four representative fruits were used for image analysis. Fruits that were positioned in the middle of the tree, at similar exposure to light, same developmental stage, and similar diameter were picked preferentially to limit intra- and inter-tree bias.

## 2.2 Estimating pigmentation pattern in red-flesh apple

### 2.2.1 Image analysis

Image acquisition, segmentation and analysis was performed as described in (Bouillon et al., 2024a). Concomitantly, a histogram-based color analysis pipeline was written in Python (Figure 1). Color analysis was inspired from a Matlab pipeline developed to study color patterning in foliar colexus (Li et al., 2021) and ColourQuant protocol (Li et al., 2022). RGB images of individual transversal apple sections were converted to gray images in order to establish three circular regions of interest (r.o.i). To be insensitive to size fruit variations, the three r.o.i of the fruit were defined as the ‘inner part’ (30-55% of the fruit radius), the ‘medium part’ (55-75%) and the outer part (75-100%). We excluded carpels (0-30%) to avoid bias in color estimation of fleshy tissue. Radius of fruits was evaluated according to the detection of centroid and edge of each fruit, from the ‘measure.regionprops’ module from scikit-image (Van der Walt et al., 2014). These circular ring masks were converted from RGB to CIEL*a*b* with Scikit-image rgb2lab function. CIEL*a*b* is a device-independent continuous color space that consists of three channels L* = “lightness”, a* = “green to magenta”, b* = “blue to yellow”. CIEL*a*b* color space enables extraction of quantitative color distribution information (Li et al., 2021). Color variations in apple fruits are constrained from green to red, so we chose the a* channel as an indicator of flesh color variations (Bouillon et al., 2024a). Normalized histograms with 10 fixed equal intervals (bins) between 0 and 60 on the a* values were then extracted for the 3 circular rings.
<figure><img src="image_5.png" /><figcaption>FIGURE 1</figcaption></figure>

Overview of the histogram-based color analysis pipeline. RGB images of individual transversal apple sections were converted to gray images in order to establish three circular regions of interest (r.o.i). To be insensitive to size fruit variations, the three r.o.i of the fruit were defined as the 'inner part' (30-55% of the fruit radius), the 'medium part' (55-75%) and the outer part (75-100%). Radius of fruit was evaluated according to the detection of centroid and edge of each fruit. These circular ring masks were converted from RGB to CIEL*a*b* color space. From a* channel, normalized histograms with 10 fixed equal intervals (bins) between 0 and 60 on the a* values were then extracted for the 3 circular rings.

## 2.2.2 Phenotypic data

Initial included 4933 images x 30 bins (10 bins for 3 circular rings) representing 492 individuals for three years (2021-2022-2023). Principal component analysis was performed on the color dataset. Correlations between coordinates of PCs and color descriptors were calculated to identify components associated with color heterogeneity/pigmentation pattern (Bouillon et al., 2024a). Each genotype was expressed as mean of its 4 fruits coordinates each year. A previous QTL analysis identified a QTL on LG16 linked to mean of a* values that colocalized with a QTL identified in a red-flesh biparental population (Bouillon et al., 2024b). In order to exclude low-colored genotypes that could display positive alleles for pigmentation pattern, thus rendering phenotyping very challenging, a subset of genotypes were filtered for mean of a* values (> 10) and the presence of LG16-QTL “F” allele at haploblock HB-16-5 (approximately 8 cM). Following this selection, 256 F1 red-flesh genotypes remained for QTL detection (Supplementary Table S1).

# 2.3 QTL analysis

## 2.3.1 Statistical analysis

Best linear unbiased prediction (BLUP) was used to estimate random effects of a mixed linear model (Piepho et al., 2007) with the following formula:

$$y = X\beta + Zu + e$$

where y is the vector of observations (phenotypic scores), β and u are vectors of fixed (year effect) and random effects (genotypic values), respectively, X and Z are the associated design matrices, and e is a random residual vector. The random effects are estimated by BLUP (Piepho et al., 2007). Because of our experimental design, not all environmental effects could be considered (tree and genotype effects were confounded), and BLUPs of across-year phenotypic values per genotype approximated the true genotypic values. BLUP was estimated from scores of individuals on the second PCA component for 2021-2022-2023 image data and used for QTL detection.

Broad-sense heritability was estimated by intra-class correlation analysis (Xu, 2022) with the following formula:

$$h^2 = \sigma_g^2/(\sigma_g^2 + \sigma_r^2)$$ where $\sigma_g^2$ and $\sigma_r^2$ were the individual genetic and residual variances, respectively.

## 2.3.2 Genotypic data

Hybrids, their parents and ancestors were genotyped using the Illumina Infinium 20K SNP array for Apple (Bianco et al., 2014). The raw SNP data was initially processed into the GenomeStudio software v2.0 and SNP data curation was performed as described by Vanderzande et al (Vanderzande et al., 2019). SNP physical position on the Malus genome (Daccord et al., 2017) was derived from the integrated genetic linkage (iGL) map for Illumina Infinium 20K SNP array (Howard et al., 2021). A total of 4220 informative SNP
were retained. SNP markers were well distributed over the 17 chromosomes. Founders were added in our pedigree from public genotyping data on GDR (Jung et al., 2018). Thus, we were able to trace marker segregation from our genotypes to the wild red-flesh accession *Malus Sieversii f. niedzwetzkyana*.

### 2.3.3 Haplotype determination

‘PediHaplotyper’ v.1.0 was used to construct haplotypes (Voorrips et al., 2016). Maximum size of each haploblock was limited to 1 cM. Identical haplotypes were described as IBD if they could be traced back to a common ancestor, whereas identical haplotypes that could not be traced back to a common ancestor were considered identical-by-state (IBS). Missing SNP data was imputed by examining both progenies and ancestors (Powell et al., 2022). Haplotypes that could not be resolved were excluded from further analysis. Finally, a total of 1104 haploblocks were considered for QTL detection.

### 2.3.4 QTL detection

QTL analyses were perfomed using the FlexQTL™ software, which implements PBA-QTL analysis via Markov Chain Monte Carlo (MCMC) simulation (www.flexqtl.nl) (Bink et al., 2014). Parameters for the analyses are reported in Supplementary Table S2. Analyses reached adequate MCMC convergence. The model applied using FlexQTL™ software was an additive, biallelic (Q/q) model where high phenotypic values were associated with Q and low phenotypic values were associated with q. Two replicate runs with different starting seed numbers were conducted to ensure repeatability of QTL results (Kostick and Luby, 2022; Powell et al., 2022). The Bayes Factor (BF; 2lnBF10) and posterior intensity values were used to examine significance and stability of a putative QTL with evidence of a QTL being considered positive, strong, or decisive if BF were 2 to < 5, 5 to < 10, and > 10, respectively (Kostick and Luby, 2022). QTL intervals were defined in a series of successive 2 cM “bins”, chromosomal segments that had a BF above 5 and the mode within a given QTL region were considered the most probable QTL position. The proportion of phenotypic variance explained (PVE) by a QTL was estimated by dividing the reported variance explained by the whole phenotypic variance (Verma et al., 2019).

Positive/negative alleles were defined via analyses of SNP haplotypes within QTL interval. Comparison of diplotypes (diploid combinations of two haplotypes) were used to infer haplotype effects. To determine if BLUP values were significantly different for presence versus absence of a given QTL haplotype, and were significantly different between the functional diplotypes, one-way analysis of variance (ANOVA) and Tukey's honest significant difference (HSD) were calculated (Powell et al., 2022). To investigate candidate gene colocalization in QTL regions, genes and TFs involved in anthocyanin synthesis were searched among QTL region using the ‘GDH13.v1.1’ apple genome (Daccord et al., 2017) as reference. Genetic positions were

defined as the boundaries of the two inner haploblocks surrounding the QTL region.

## 2.4 Comparative analysis between RF and WF

### 2.4.1 Plant material

To compare RF and WF, a red-flesh parent from the population, cv. ‘R201’ (plantation year: 2013) that exhibited pigmentation pattern was chosen. Red-flesh fruit are fully red after fruit set (Espley et al., 2012). Fruit flesh pigmentation pattern develop during fruit development and is stable through maturity. To identify the biochemical and genetic factors involved in anthocyanin depigmentation, the transition stage between fully red fruits and heterogeneously pigmented fruits was targeted during fruit development (84 Days After Full Bloom - DAFB). pH was measured on red and white-flesh fruit sections using a pH-meter Seven Compact S210 (Metler Toledo) as described in (Bouillon et al., 2024a).

### 2.4.2 Determination of phenolic contents in RF and WF

Methanol (CHROMASOLV for LC-MS), acetic acid, formic acid and dimethylsulfoxide (analytical reagent grade) were purchased from Fisher Scientific. Ultrapure water was obtained from a MilliQ advantage A10 purification system (Millipore). Hyperoside, cyanidin-3-O-galactoside chloride, cyanidin-3-O-arabinoside chloride, procyanidin B1, procyanidin B2, procyanidin C1, (+)-catechin, (-)-epicatechin, chlorogenic acid, phlorizin, avicularin and quercitrin standards were purchased from Extrasynthese. 4-p-coumaroylquinic acid was purchased from Ambinter. Reynoutrin was purchased from Clinisciences. Daidzein (internal standard) was purchased from Molekula. All standards were furnished with a certificate of analysis.

Stock solutions of procyanidin B1, procyanidin B2, procyanidin C1, phlorizin, hyperoside, avicularin, quercitrin, reynoutrin and daidzein were prepared in DMSO at a concentration of 5 mg.mL⁻¹. Stock solutions of (+)-catechin and (-)-epicatechin were prepared in MeOH at a concentration of 5 mg.mL⁻¹. Stock solutions of chlorogenic acid and 4-p-coumaroylquinic acid were prepared in H₂O at a concentration of 5 mg.mL⁻¹. Stock solution of cyanidin-3-O-galactoside and cyanidin-3-O-arabinoside chloride were prepared in MeOH containing 5% acetic acid at a concentration of 5 mg.mL⁻¹. All the dilutions were then carried out in MeOH except for cyanidin-3-O-galactoside and cyanidin-3-O-arabinoside chloride in MeOH containing 5% acetic acid. Each calibration curves were set to the concentration range expected for each compound. The linearity of the calibration curves was assessed by injecting 6 levels of calibration standards in three replicates. Residuals and their distribution were monitored. LOD and LOQ were defined as 3- and 10- fold the signal-to-noise, respectively.
Precision was evaluated by repeated analysis of standards within different analytical batches. Each sample were spiked with the same level of daidzein at the start of extraction procedure in order to evaluate the recovery of extraction. Data acquisition and analysis were conducted on Microsoft Excel 365® and the Minitab 19® software.

Fruits without peel were cuted and two flesh pieces were selected on the sun-exposed side of the fruit (about 5 grams per fruit). Samples were immediately frozen in liquid nitrogen, stored at -80°C and freeze-dried. This dried material was crushed into a fine and homogeneous powder and stored until analysis. Polyphenols extraction protocol used was previously described by (Bouillon et al., 2024a).

Quantitative LC-UV analysis was performed on a ThermoFisher Scientific Vanquish Flex UHPLC system equipped with a quaternary solvent manager, a sample manager and a column compartment as chromatographic part and coupled to a variable wavelength detector as detection part. A sample volume of 10 µL was injected onto a ZORBAX RRHD StableBond Aq column (2.1 x 150 mm; 80 Å; 1.7 µm). Samples were kept at 10°C and column was maintained at 30°C with a flow rate of 0.3 mL.min-1. The mobile phase consisted of 0.1% formic acid in both H2O (A) and MeOH (B) used in gradient mode as follows: from 5 to 15% B in 5 min, then 15 to 35% B from 5 to 30 min, then 35 to 50% B from 30 to 35 min, then 50 to 100% B from 35 to 37 min, hold at 100% B from 37 to 39 min, and afterward the column was re-equilibrated at initial conditions during 8 min. The UV detection was measured at different wavelengths optimized for individual compounds and

are listed in Table 1. Data were processed using Chromeleon 7.3.2 software.

### 2.4.3 RNAseq analysis and GO-enrichment analysis

Six samples of 'R201' with 3 semi-biological replicates (fruits from different trees) and distinction of cortical WF and RF were analyzed by RNA sequencing. Total RNA was extracted from fruit flesh according to the method detailed by Celton et al (Celton et al., 2014). The library construction and Pair Ends 100bp reads sequencing were completed by BGI using Illumina technology, and the quality was checked using FastQC (Andrews, 2010). Alignement, quantification and differential expression analysis were performed with 'AnaDiff' (Pelletier, 2022). Lists of DEGs were obtained with 'DESeq2' (Love et al., 2014) using a p-value and BH, a method for controlling the false discovery rate (FDR), of 0.05. The control represents the WF section while the treatment is the RF outer flesh. For graphical display, colorblindness-friendly colors were provided by 'viridis' package.

A new functional annotation of GDDH13v1.1 (Daccord et al., 2017) was performed according to the latest proteomic Arabidposis thaliana homology. The blast tool was used with a 0.05 p-value threshold. *Arabidopsis thaliana* information were extracted from TAIR10 database (Cheng et al., 2017). Interproscan (Quevillon et al., 2005) was also used on GDDH13 v1.1 proteome. A compilation of blast results, TAIR10 database and Interproscan results enabled obtention of the new functional annotation of *Malus domestica*. From the RNAseq results, we selected differentially

<table><caption>TABLE 1 Summary of validation parameters.</caption><thead><tr><th>Compound</th><th>Retention time (min)<sup>a</sup></th><th>UV detection max (nm)</th><th>LOD (ng)<sup>b</sup></th><th>LOQ (ng)<sup>c</sup></th><th>Calibration range (ng)<sup>d</sup></th><th>R<sup>2</sup></th></tr></thead><tbody><tr><td>Procyanidin B1</td><td>12.75 ± 0.12</td><td>279</td><td>0.22</td><td>0.44</td><td>0.9 - 44.1</td><td>0.992</td></tr><tr><td>(+)-catechin</td><td>13.45 ± 0.34</td><td>279</td><td>0.5</td><td>1</td><td>2.5 - 99.8</td><td>0.999</td></tr><tr><td>Procyanidin B2</td><td>16.26 ± 0.07</td><td>279</td><td>0.25</td><td>0.5</td><td>9.9 - 496.5</td><td>1</td></tr><tr><td>(-)-epicatechin</td><td>16.74 ± 0.25</td><td>279</td><td>0.25</td><td>0.5</td><td>5 - 248</td><td>0.999</td></tr><tr><td>Chlorogenic acid</td><td>17.44 ± 0.39</td><td>325</td><td>0.49</td><td>0.99</td><td>494.5 - 4945</td><td>0.999</td></tr><tr><td>Cyanidin-3-O-galactoside chloride</td><td>18.62 ± 0.11</td><td>512</td><td>2.4</td><td>4.7</td><td>4.7 - 236.3</td><td>0.999</td></tr><tr><td>4-p-coumaroylquinic acid</td><td>20.79 ± 0.34</td><td>312</td><td>0.13</td><td>0.25</td><td>5 - 250</td><td>0.999</td></tr><tr><td>Procyanidin C1</td><td>21.99 ± 0.33</td><td>279</td><td>0.87</td><td>2.2</td><td>4.3 - 216.8</td><td>0.999</td></tr><tr><td>Cyanidin-3-O-arabinoside chloride</td><td>23.95 ± 0.12</td><td>512</td><td>2.4</td><td>4.7</td><td>4.7 - 235.3</td><td>0.998</td></tr><tr><td>Phloridzin</td><td>33.41 ± 0.11</td><td>284</td><td>0.09</td><td>0.23</td><td>23.1 - 692.3</td><td>0.998</td></tr><tr><td>Quercetin 3-O-galactoside</td><td>36.21 ± 0.17</td><td>350</td><td>0.24</td><td>0.47</td><td>0.5 - 9.5</td><td>0.999</td></tr><tr><td>Quercetin 3-O-xyloside</td><td>36.78 ± 0.04</td><td>350</td><td>0.1</td><td>0.24</td><td>0.2 - 7.3</td><td>0.998</td></tr><tr><td>Quercetin 3-O-arabinoside</td><td>37.33 ± 0.03</td><td>350</td><td>0.09</td><td>0.24</td><td>0.2 - 7.1</td><td>0.995</td></tr><tr><td>Quercetin 3-O-rhamnoside</td><td>37.65 ± 0.01</td><td>350</td><td>0.23</td><td>0.47</td><td>2.3 - 93.2</td><td>0.999</td></tr><tr><td>Daidzein (IS)</td><td>38.73 ± 0.19</td><td>250</td><td>0.1</td><td>0.25</td><td>2.5 - 49</td><td>0.996</td></tr></tbody></table>

<footnote><sup>a</sup>Values are mean ± SD (n=3); <sup>b</sup> The signal-to-noise was set to 3:1; <sup>c</sup>The signal-to-noise was set to 10:1; <sup>d</sup>Suitable for samples.</footnote>
expressed genes with a BH < 0.05 and logFC > |0.5|. If logFC > 0.5, genes were considered as 'RF overexpressed' in comparison to WF. On the contrary, if the logFC < -0.5, genes were considered as 'RF underexpressed' in comparison to WF. GO enrichment was computed for 'RF overexpressed' and 'RF underexpressed' genes. These two lists of genes were treated identically but distinctly for a GO enrichment by an hypergeometric law (0.05 p-value threshold) against all the genes of the genome. GO visualizations for the three GO types (Biological Process, Molecular Function and Cellular Component) were performed with Plotly (Inc., 2015).

### 2.4.4 Differentially methylated regions (DMRs) analysis

The same six 'R201' samples used for differential RNAseq analyses were used. DNA was purified according to the CTAB protocol detailed by (Porebski et al., 1997). Library construction and whole genome bisulphite sequencing (PE 100 bp) were performed by BGI TECH SOLUTIONS (Hong Kong). The quality of paired end sequencing reads (100 bp) generated from WGBS was evaluated with FastQC (version 0.11.9). The Trimmomatic (v0.39) software was used to remove adaptor sequences and poor quality bases using the given parameters "LEADING:3, TRAILING:3, MINLEN:80, ILLUMINACLIP:adapter:2:30:10". The reference genome GDDH13 v1.1 was converted (C to T and G to A) by 'bismark_genome_preparation' tool within Bismark (v0.24.0) (Krueger and Andrews, 2011). Subsequently, trimmed sequencing reads were aligned to the converted genome version using bowtie2. The resulting aligned BAM files were then processed to remove duplicates with the deduplicate_bismark tool. Methylation extraction was performed on the deduplicated reads with bismark_methylation_extractor tool to generate the context-specific methylation-states (CpG, CHG and CHH) at the individual cytosine-level for each sample with the following parameters "-p -parallel 4 -bedGraph -CX_context -cytosine_report -comprehensive -gzip". From this extraction, a genome-wide cytosine report was generated.

The R package methylkit (v1.26.0) (Akalin et al., 2012) was used to compare the methylomes between red and white fruit flesh. Pairwise comparison was conducted with two and three semi-biological replicates, respectively. Methylation levels at cytosines in CpG, CHG or CHH contexts were initially calculated using the 'methRead' function with mincov=4. The 'tileMethylCounts' function was used to summarize methylated and unmethylated base counts over tiling windows across the genome with the parameters 'win.size=100, step.size=100, cov.bases=4'. The 'unite' function was used to merge and retain only bases that were covered in tiled regions of all replicates with the destrand = FALSE parameter. Subsequently, the differential methylation statistics were calculated with 'calculateDiffMeth' function and 'getMethylDiff' function was used to select the DMRs for CpG, CHG, and CHH contexts that satisfied the following parameters 'minDiff=25 q-value=0.05'. The command 'bedtools closest' from Bedtools (v2.27.1) (Quinlan and Hall, 2010) with the parameter '-Db' was employed on the GDDH13 annotated genes to identify and retrieve the nearest gene to each DMR (Supplementary Table S7).

## 3 Results

### 3.1 Principal component analysis and color descriptors relationships

Correlation analysis between color descriptors (Bouillon et al., 2024a) and principal components (PCs) confirmed correlation between first axis and mean of a* values (Figure 2A) with a Pearson correlation value of 0.8. Second axis was positively correlated to standard deviation of a* (Pearson correlation = 0.5). As expected, PC2 captured variance associated with red-flesh spatial distribution and was considered as the phenotypic trait defining pigmentation pattern (Figure 2B). PC1 and PC2 were not correlated, allowing decorrelation of red-flesh intensity and distribution.

Broad-sense heritability of pigmentation pattern was 0.56 and confirmed the potential genetic control of this trait. BLUP pigmentation pattern varied from -4.446 to 5.111 (Figure 2C) and achieved distinction between homogeneous pigmentation (white-off/yellow and pink/red homogeneous pigmented phenotypes are confounded) and pigmentation pattern (Figure 2D).

### 3.2 QTL detection and haplotype analysis

To identify genetic regions that significantly affect pigmentation pattern, PBA-QTL analysis was performed with BLUP of pigmentation pattern used as phenotypic trait. One QTL for BLUP pigmentation pattern was detected on LG8 (Figure 3A; Supplementary Table S3). This QTL region was located between 23 and 42 cM with QTL mode located at 30 cM. Phenotypic variance explained (PVE) by this QTL was 20.2% (Supplementary Table S4). Haplotype analysis was conducted on LG8-QTL (Figure 4). Haploblock was constructed from six SNPs (Figure 3B) between 30.04 and 30.67 cM. Five haplotypes were identified in this population and noted H11, H12, H13, H19 and H24 (Figure 3B). Combination of diplotypes exhibited significant phenotypic variations on pigmentation pattern (p < 0.05; Figure 4A). Genotypes that displayed at least one H11 haplotype showed significantly higher phenotypic mean for pigmentation pattern (p < 0.05; Figure 4B). The red-flesh parent 'R201' was heterozygous (H19/H11) at this QTL and displayed one positive haplotype associated with pigmentation pattern.

Candidate genes associated with the LG8-QTL region were investigated (Figure 3A). Four putative candidate genes (among 724 genes) were identified within that region: two genes annotated for MYB (MD08G1092000 and MD08G1117700), WRKY11 (MD08G1127200) and HY5 (MD08G1147100) were located at 23, 30.6, 32.48 and 37.97 cM respectively.

### 3.3 Wide-targeted metabolite analysis between red- and white-flesh

To determine if there was significant differences in phenolic compounds accumulation between RF and WF, wide-targeted
<figure><img src="image_6.png" /><figcaption>FIGURE 2</figcaption></figure>

Principal Component Analysis (PCA) for the Pigmentation Patterning defined by color density based analysis. (A) Pearson correlation between principal components and color descriptors. (B) PCA plot. Genotypes are colored according to the standard deviation of a* values. (C) Histogram of values for the 2nd principal component (phenotypic trait of pigmentation patterning). (D) Visual examples of values for the 2nd principal components for different ranges of a* means.

metabolite data were compared for each main phenolic groups: dihydrochalcones (phloridzin), hydroxycinnamic acids (4-p-coumaroyl quinic acid, chlorogenic acid), flavonols (quercetin 3-O-galactoside, quercetin 3-O-xyloside, quercetin 3-O-arabinoside and quercetin 3-O-rhamnoside), flavan-3-ols (catechin, epicatechin, procyanidin B1, procyanidin B2 and procyanidin C1)

and anthocyanins (cyanidin 3-O-galactoside and cyanidin 3-O-pentoside).

The minor anthocyanin of the fruit was quantified as equivalent of cyanidin-3-O-arabinoside (Table 1). Indeed, an UHPLC-HRMS analysis of the extract confirmed the molecular formula C20H19O10, indicating the presence of a cyanidin pentoside.

<figure><img src="image_7.png" /><figcaption>FIGURE 3</figcaption></figure>

Bayes Factor (BF) across the LG8 and genotypes for LG8-QTL haploblock (A) BF across the LG8 showing significant association with pigmentation pattern. The horizontal red dashed line indicated positive threshold. The putative candidate genes underpinning the QTL region are shown. (B) The five haplotypes associated with LG8-QTL haploblock (HB-08-26) and their SNP sequences.

B

<table><caption>LG8 - QTL mode (SNPs index)</caption><thead><tr><th>Haplotype</th><th>F1 frequency</th><th>2093</th><th>13349</th><th>13351</th><th>13352</th><th>13353</th><th>13354</th></tr></thead><tbody><tr><td>Genetic position (cM)</td><td></td><td>30.044</td><td>30.09</td><td>30.6</td><td>30.6</td><td>30.61</td><td>30.61</td></tr><tr><td>Physical position (GDDH13v1v.1)</td><td>10122497</td><td>10140125</td><td>10612716</td><td>10613040</td><td>10615421</td><td>10617538</td><td></td></tr><tr><td>H11</td><td>340</td><td>G</td><td>T</td><td>G</td><td>C</td><td>C</td><td>A</td></tr><tr><td>H12</td><td>174</td><td>G</td><td>C</td><td>G</td><td>C</td><td>C</td><td>A</td></tr><tr><td>H13</td><td>169</td><td>G</td><td>C</td><td>A</td><td>T</td><td>T</td><td>G</td></tr><tr><td>H19</td><td>116</td><td>G</td><td>T</td><td>A</td><td>T</td><td>T</td><td>G</td></tr><tr><td>H24</td><td>150</td><td>A</td><td>T</td><td>G</td><td>C</td><td>C</td><td>A</td></tr></tbody></table>
<figure><img src="image_8.png" /><figcaption>FIGURE 4</figcaption></figure>

Distribution of BLUP pigmentation pattern for the five F1 offsprings haplotypes for LG8-QTL. (A) Distribution of the 11 diplotypes for LG8-QTL. (B) Distribution of the five F1 offsprings according to the absence/presence of H11 haplotype. The median (denoted by a horizontal bar in the box), the 25th percentile (denoted by the bottom edge of the box), the 75th percentile (denoted by the top edge of the box), the mean (denoted by a black dot) and the dots indicate single observations. Number of genotypes are listed above. Significantly different phenotypic means between segregating classes are identified by different letters (tukey HSD, p < 0.05).

Three different cyanidin pentosides have been reported in apple: cyanidin-3-O-arabinoside, cyanidin-3-O-xyloside and cyanidin-7-O-arabinoside, the first one being the most cited (Slatnar et al., 2014; Kokalj et al., 2019). Unfortunately, the retention time does not match and the quantification in equivalent was chosen given that cyanidin-3-O-arabinoside standard was co-injected. WF accumulated significantly more dihydrochalcones and hydroxycinnamic acids with an average concentration of 495.1 and 20755.6 µg/g of dry weight (DW), respectively, against 287.6 and 9603.8 µg/g of DW in RF (Figure 5; Supplementary Table S8). Data on flavonols were more balanced, nonetheless, RF contents were higher (64.2 µg/g of DW) in comparison to WF (52.6 µg/g of DW). RF accumulated significantly more flavan-3-ols and anthocyanins with an average concentration of 761.1 and 484.7 µg/g of DW, respectively, compared to 586.8 and 37.7 µg/g of DW in WF (Figure 5).

Additionally, pH measurements were also carried out on RF and WF. RF showed lower pH values with a mean of 2.97 against 3.11 for WF (Supplementary Table S5).

## 3.4 Analysis of differentially expressed genes between RF and WF by RNA-seq and GO-functional analyses

To functionally characterize flesh pigmentation pattern, we performed metabolomic and transcriptomic comparisons between RF and WF (Figures 5, 6). The transition step between fully red fruits and heterogeneously pigmented fruit were targeted during fruit development (84 Days After Full Bloom - DAFB). A total of 2,727 DEGs were identified using a false discovery rate (BH) < 0.05. A total of 853 genes were up-regulated (log Fold Change logFC > 0.5) in RF and 718 genes were down-regulated (logFC < -0.5) in RF compared to WF. Among the anthocyanin-related genes, MdGSTU22 (MD17G1272100), MdCHS (MD04G1003300),

MdPAL (MD04G1096200) and MdDFR (MD15G1024100) had the highest differential expression between RF and WF sections (Figure 6). This analysis included three semi-biological replicates of RF and WF (fruits from different trees) that are fully reproducible (Supplementary Figure S1).

To explore the biological processes of DEGs, we classified the gene ontology (GO) function and analyzed the enrichment (Supplementary Figure S2). GO terms were divided into three functional categories: molecular function, cellular component and biological process. DEGs up-regulated in RF were enriched in GO terms associated with response to environmental stimuli ('response to light stimulus', 'response to cold', 'response to temperature stimulus', 'response to UV') and flavonoid biosynthetic process ('regulation of anthocyanin biosynthetic process', 'flavonoid biosynthetic process', 'anthocyanin-containing compound metabolic process'). DEGs up-regulated in WF were enriched in GO terms associated with hormonal signal transduction ('response to gibberellin', 'response to brassinosteroid') and flavonoid biosynthetic process ('flavonoid biosynthetic process', 'flavone synthase activity').

## 3.5 Expression of anthocyanin biosynthetic genes

To further investigate transcriptional differences between RF and WF, DEGs encoding TFs related to MBW complex, early anthocyanin biosynthesis gene (EABG) and late anthocyanin biosynthesis gene (LABG) were analyzed. A total of 32 DEGs related to anthocyanin biosynthesis were up-regulated in RF compared to WF with logFC > 0.5 (Figure 7).

Globally, genes related to anthocyanin biosynthesis were up-regulated while genes related to phenolic compounds and anthocyanin catabolism were down-regulated in RF compared to WF (Figure 7).
<figure><img src="image_9.png" /><figcaption>FIGURE 5</figcaption></figure>

Concentration of major phenolic compounds between red-flesh (RF) and white-flesh (WF). The concentrations of these metabolites are expressed in µg.g⁻¹ of dry weight. Phenolic compounds biosynthesis (enzymes and precursors) are represented. PAL, phenylalanine ammonia lyase; C4H, cinnamate 4-hydroxylase; 4CL, 4-coumarate coenzyme A ligase; CHS, chalcone synthase; CHI, chalcone isomerase; F3H, flavanone 3-hydroxylase; F3'H: flavonoid 3'-hydroxylase; DFR, dihydroflavonol 4-reductase; ANS, anthocyanidin synthase; UFGT, (UDP)-glucose:flavonoid 3-O-glucosyltransferase; FLS, flavonol synthase; LAR, leucoanthocyanidin reductase; ANR, anthocyanidin reductase. The mean (black dot) and the standard deviation (vertical black line) are represented in each range-plot for RF and WF. Significantly different means are identified by different letters (tukey HSD, p < 0.05).

## 3.6 Differential methylation analysis
Between RF and WF

To explore a potential cause of the observed differential RF-WF gene expression, DNA methylation analysis was performed on the same tissues as those used for the RNAseq analysis. DNA methylation may have an important regulatory effect on gene expression (Zhang et al., 2018). Methylome analysis of RF and WF (Supplementary Figure S3) permitted the identification of 1,320 differentially methylated regions (DMRs) in the three contexts (CG,

CHG, and CHH), comprising 270 hypomethylated and 1,050 hypermethylated regions in RF compared to WF. Concomitantly, 73 DMRs were associated with DEGs and 2 were found in the proximity of genes involved in anthocyanin synthesis: *MdGSTU22* (MD17G1272100) and *MdPAL1* (MD04G1096200), with a difference in the percentage of methylation of -32.3% and -34.8% respectively. The DMR (CHH context) associated to *MdGSTU22* was located 1,862 bp downstream of the 3' UTR of the gene. The DMR (CHG context) associated to *MdPAL1* was located at 1,514 bp upstream of the 5' UTR of the gene.
<figure><img src="image_10.png" /><figcaption>FIGURE 6</figcaption></figure>

Volcano plot of differentially expressed genes (DEGs) between RF and WF The negative log of p-value (base 10) is plotted on the Y-axis, and the log of the FC (base 2) is plotted on the X-axis. The purple points on this graph represent genes that are significantly differentially expressed between red-flesh (RF) and white-flesh (WF). A total of 2,727 DEGs were identified using a false discovery rate (BH) < 0.05. A total of 853 genes were up-regulated (log Fold Change logFC > 0.5) in RF and 718 genes were down-regulated (logFC < -0.5) in RF compared to WF. The gray points represent genes that are not differentially expressed with log FC is > |0.5| and BH < 0.05.

# 4 Discussion

## 4.1 Genetic control of fruit flesh pigmentation pattern

The genetic determinism of anthocyanin biosynthesis in apple is well characterized and is mainly regulated by the MBW complex through *MdMYB10* activity (Honda and Moriya, 2018; Wang and Chen, 2021). Recently, new regulators of anthocyanin biosynthesis have been identified (Liu et al., 2021), but less is known about spatial regulation of anthocyanin deposition in the fruit flesh. Here we have shown that color patterning of apple red flesh is under genetic control. Indeed, fruit flesh pigmentation patterns were found heritable over three years (2021-2023) with a heritability value of 0.56, thus enabling QTL detection. A QTL for BLUP pigmentation pattern was identified on LG08 and colocalized with two genes that encoded MYB TFs. One R2R3-MYB TF is localized upstream of the QTL mode in the LG8-QTL region and annotated as *MdMYB44* (approximately 23 cM). One uncharacterized MYB-3R colocalized with the LG8-QTL mode (approximately 30.6 cM).

It has been proposed that R2R3-MYB evolved from MYB-3R after the loss of the R1 repeat sequence (Feng et al., 2017). MYB dual activity are well known in anthocyanin synthesis repression (LaFountain and Yuan, 2021) and some R2R3-MYB anthocyanin repressors have been identified in *Malus domestica* (Wang and Chen, 2021). MYB repressors may inhibit the MBW complex via the competition with MYB activators for bHLH binding sites (Yan et al., 2021). Interestingly, *MdMYB44* has been linked to color fading in crabapple leaves by downregulating ABGs (Meng et al., 2022). Recent work suggests that *MdMYB44* negatively regulates vacuolar acidification and decreases malate content (Jia et al., 2021). Vacuolar pH is determinant in anthocyanin color and stability (Zhang et al., 2014). pH variations within apple flesh could lead to anthocyanin-related color loss and differential patterning. Allelic variations in R2R3-MYB sequences are often associated with

important phenotypic variations (Espley et al., 2009; Castillejo et al., 2020; Jia et al., 2021) and allelic variations in a new MYB anthocyanin-related activator or repressor could determine formation of pigmentation pattern in red-flesh apple.

Genetic factors other than the genes involved in the MBW complex can regulate anthocyanin biosynthesis (Wang and Chen, 2021). *MdWRKY11* increased the accumulation of flavonoids and anthocyanin in RF apple by promoting LABGs (F3H, FLS, DFR, ANS, and UFGT genes) expression (Wang et al., 2018). However, *MdWRKY11* does not interact with the MBW complex but interacts with *MdHY5* and promotes its activity (Liu et al., 2019). Indeed, the photoresponse factor *MdHY5* was shown to promote anthocyanin accumulation by regulating the expression of the *MdMYB10* gene and downstream ABGs (An et al., 2017). As a consequence, allelic variations associated with *MdWRKY11* and *MdHY5* could be involved in flesh pigmentation patterning.

## 4.2 Developmental cues in pattern acquisition

Whether red-flesh fruits are homogeneously or heterogeneously pigmented at harvest, they are fully red after fruit set and pattern formation only occurs during fruit development. The transition step between fully pigmented fruits and heterogeneously pigmented fruit is determinant in spatial anthocyanin deposition. Wide-targeted metabolite analysis identified differential accumulation in phenolic compounds between RF and WF. As expected, anthocyanins were 15-fold more accumulated in RF than in WF. This drastic difference coincided with increased gene expression of ABGs and MBW complex. Interestingly, some ohnologous genes (Lallemand et al., 2023) exhibited similar differential expression patterns between RF and WF.

The regulatory genes *MdMYB10* and *MdbHLH3* were up-regulated in RF compared to WF (Figure 7). MYB-TFs
<figure><img src="image_11.png" /><figcaption>FIGURE 7</figcaption></figure>

Expression pattern of genes involved in anthocyanin biosynthesis, degradation, and transport pathway between red-flesh and white-flesh. Color boxes from left to right represent values of log Fold Change. PAL, phenylalanine ammonia lyase; 4CL, 4-coumarate coenzyme A ligase; C4H, cinnamate 4-hydroxylase; CHS, chalcone synthase; CHI, chalcone isomerase; F3H, flavanone 3-hydroxylase; F3'H: flavonoid 3'-hydroxylase; DFR, dihydroflavonol 4-reductase; ANS, anthocyanidin synthase; UFGT, (UDP)-glucose:flavonoid 3-O-glucosyltransferase; FLS, flavonol synthase; LAR, leucoanthocyanidin reductase; ANR, anthocyanidin reductase; GST, glutathione S-transferase; ABC, ATP-binding cassette; LAC, laccase; PPO, polyphenol oxidases; POD, class III peroxidases. Color gradient, ranging from yellow to purple represents low, middle, and high values of differential gene expression.

anthocyanin activators are able to bind directly to the promoter of the key gene for anthocyanin biosynthesis and therefore, affect the final anthocyanin accumulation (Yan et al., 2021). As a consequence, positive regulation by MYB is generally associated with a global up-regulation of ABGs (Espley et al., 2012). An additional level of anthocyanin regulation involves WD40 proteins, formerly identified as TTG1 in A. thaliana (Gonzalez et al., 2007; Brueggemann et al., 2010), which interact with bHLH through the MBW complex (An et al., 2012). MdTTG1 (MD01G1228700) was expressed through WF and RF but did not show differential expression (logFC = -0.26). This suggests a required interplay between WD40, MYB-TF and bHLH to upregulate anthocyanin synthesis. In addition to MdMYB10, other MYB-TFs were differentially expressed between WF and RF. MYB3 is a transcriptional repressor involved in repression of phenylpropanoid and anthocyanin biosynthesis in A. thaliana (Zhou et al., 2017; Kim

et al., 2022). DEGs analysis identified a MdMYB3 (MD15G1411200), down-regulated in RF (logFC = -1.9), making it a potential anthocyanin repressor in red-flesh apple. MdMYB9 is a positive regulator of anthocyanin and proanthocyanidin biosynthesis in red-flesh apple through the HY5-miR858-MYB9 loop (Li et al., 2023) and jasmonate-mediated signal (An et al., 2021). Subsequent functional analyses demonstrated that MdMYB9 bound to the promoter regions of flavonoid biosynthetic genes and induced anthocyanin accumulation in calli in the presence of JA (Honda and Moriya, 2018). Unexpectedly, MdMYB9 was down-regulated in RF (logFC = -0.7, BH = 0.007) while MdHY5 was not differentially expressed between WF and RF, challenging the current MYB9 related regulatory mechanism. Nonetheless, DEGs related to jasmonate were over-represented in RF, which support the pivotal role of jasmonate signaling in HY5-miR858-MYB9 loop involved in anthocyanin regulation. Plant hormones can regulate anthocyanin biosynthesis
through different pathways (Shi et al., 2023). GO-enrichment analysis showed that 'response to auxin' related DEGs were over-represented in RF and WF. Interestingly, DEGs related to brassinosteroid and gibberelin were over-represented in WF while DEGs related to jasmonate were over-represented in RF (Supplementary Figure S2). Further studies based on hormone measurements, miRNAs and gene expression analysis during fruit development kinetics could help to identify new regulators, mediated by different signals, of anthocyanin biosynthesis. Interestingly, flavan-3-ols and flavonols were also higher in RF (Figure 5). Flavonols are derivates of dihydroflavonols synthesized by flavonol synthase (FLS) and flavan-3-ols are synthesized from leuco/anthocyanidins by LAR and ANR. MdANR (MD01G1077800, MD05G1335600), MdLAR (MD06G1211400, MD13G1046900, MD16G1048500) and MdFLS (MD08G1168600, MD15G1353800) were not differentially expressed between RF and WF (Supplementary Table S6). Inversely, dihydrochalcones and hydroxycinnamic acids were more abundant in WF than in RF (2-fold). Differential accumulation of early phenolic compounds could be associated with remobilization of flavonoid metabolism or with products resulting from anthocyanin catabolism given that enzymes related to anthocyanin catabolism such as, laccases (LAC), polyphenol oxydases (PPO) and peroxydases (POD), were more expressed in RF compared to WF. Moreover, anthocyanins' color and stability depend on vacuolar pH via their pH-dependent forms (Zhang et al., 2014). RF showed lower pH values (pH = 2.97) than WF (pH = 3.11), therefore contributing to anthocyanin degradation in WF.

The translocation of glycosylated anthocyanins into the vacuole constitutes a limiting step in anthocyanin accumulation (Buhrman et al., 2022). GSTs are involved in anthocyanin transport by conjugating to ABC transporters (Buhrman et al., 2022). A previous study identified MdGSTU12 as being involved in apple anthocyanin accumulation (Zhao et al., 2021b). Among the genes with high expression levels (Figure 6), MdGGSTU22 was up-regulated in RF (logFC = 2.2) compared to WF suggesting that this gene is related to anthocyanin accumulation in RF. A recent study (Eichenberger et al., 2023), showed that arGST plays a catalytic role in anthocyanin biosynthesis. However, we did not observe flavonoid substrate competition between RF and WF. Further functional validation is required to understand MdGSTU22 implication in anthocyanin biosynthesis.

Among the 1,320 DMRs (Supplementary Figure S3), we identified a DMR associated with MdGSTU22 (-32.3% of methylation) as observed in red-flesh kiwi-fruit (Zhang et al., 2023). This DMR (CHH context) was located 1,862 bp downstream of the 3' UTR of the gene. Promoter regions of genes are generally located in the 5' region. However, it happens that cis-regulatory elements may also be located in the 3' region and affect gene expression (Bernardes and Menossi, 2020). Since the function of those 3' regulatory region in apple are not well-known, a DMR in this region might affect the expression of this particular gene. One ABG (MdPAL1) at the start of flavonoid pathway was also linked to a DMR (-34.8% of methylation). These two hypomethylated DMRs were associated with DEGs up-regulated in RF compared to WF, suggesting that DNA methylation could modify anthocyanin-related genes expression as observed in apple skin (Telias et al.,

2011; Jiang et al., 2019; Xu et al., 2024). While pigmentation is mainly attuned to developmental process, pigment production and accumulation is also dependent on environmental stimuli (especially light and temperature) (Espley and Jaakola, 2023). GO-enrichment analysis confirmed that RF was more dependent to environmental stimuli than WF with enrichment in GO term related to light ('response to UV', 'response to light stimulus') and temperature ('response to cold', 'response to temperature stimuli'). Further investigations should be conducted to estimate the genotype, environmental stimuli and hormone related signals effects (Shi et al., 2023), and their interaction, on fruit flesh pigmentation patterning.

## 4.3 Color fading drives pigmentation heterogeneity

QTL analysis revealed for the first time a new locus associated with the genetic determinism of red-flesh pigmentation pattern in apple. The fact that we did not identify a DEG that colocalized within the QTL region could be explained by the experimental setup where we targeted one developmental stage (transition stage between fully red fruits and heterogeneously pigmented fruits) while other activator/repressor of anthocyanin biosynthesis could also be expressed before. QTL may contain master regulators that can control the expression of many genes (Shen et al., 2022). Moreover, PVE of this QTL was 20.2% which suggests that other genetic factors could be involved in pigmentation pattern. This study identified two potential mechanisms involved in pigmentation patterning in red-flesh apple: an activator-repressor system, which may be genetically inherited, as suggested by QTL detection, and color fading, which can occur during fruit development. These mechanisms are intricated, and likely influenced by developmental and environmental factors. Broader QTL analyses with progenies that segregate between homogeneous and heterogeneous pigmented phenotypes, or Genome Wide Association Studies (GWAS) could help to fine mapping LG8-QTL, and improve our understanding of the genetic determinism associated with this original phenotypic trait.

Altogether, our results support that WF result from color fading while maintenance through development of RF is associated with an increased expression of anthocyanin-related positive regulators. Indeed, while ABGs and MBW complexes were up-regulated in RF, a lack of anthocyanin transport (down-regulation of MdGSTU22 in WF compared to RF) and an increase in phenolic catabolism (up-regulation of MdLAC, MdPPO and MdPOD in WF compared to RF) were concomitant with the color fading process in WF. Color fading has been identified in apple organs such as flower (Han et al., 2020), leaves (Meng et al., 2022) and fruit skin (Xu et al., 2024). Here, color fading was linked to variations in anthocyanin related genes expression while methylation levels could modulate inner pigmentation pattern. Extension of comparative study on other genotypes could enable the identification of new anthocyanin regulatory mechanisms involved in anthocyanin spatial deposition and therefore flesh pigmentation pattern.
## Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material. GEO accession numbers (online repositories): GSE269807 and GSE269214.

## Author contributions

PB: Writing – original draft, Writing – review & editing. EB: Writing – review & editing. A-LF: Writing – review & editing. SB: Writing – review & editing. SH: Writing – review & editing. YL: Writing – review & editing. MC: Writing – review & editing. FF: Writing – review & editing. AB: Writing – review & editing. DB: Writing – review & editing. FB: Writing – review & editing. J-MC: Writing – review & editing.

## Funding

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This work was financed by ANRT, IFO and IRHS (CIFRE convention N° 2021/0182).

## Acknowledgments

We would like to thank Dr. Mao Li for advices about color patterning analysis. The authors acknowledge PTM PHYTO (SFR QUASAV) for the instruments used for determination of phenolic content. The authors acknowledge PTM ANAN for RNA and DNA quality check.

## References

Akalin, A., Kormaksson, M., Li, S., Garrett-Bakelman, F. E., Figueroa, M. E., Melnick, A., et al. (2012). methylkita: a comprehensive r package for the analysis of genome-wide dna methylation profiles. *Genome Biol.* 13, R87. doi: 10.1186/gb-2012-13-10-r87

Allan, A. C., Hellens, R. P., and Laing, W. A. (2008). Myb transcription factors that colour our fruit. *Trends Plant Sci.* 13, 99–102. doi: 10.1016/j.tplants.2007.11.012

An, J.-P., Qu, F.-J., Yao, J.-F., Wang, X.-N., You, C.-X., Wang, X.-F., et al. (2017). The bzip transcription factor mdhy5 regulates anthocyanin accumulation and nitrate assimilation in apple. *Horticulture Res.* 4. doi: 10.1038/hortres.2017.23

An, J., Xu, R., Liu, X., Zhang, J., Wang, X., You, C., et al. (2021). Jasmonate induces biosynthesis of anthocyanin and proanthocyanidin in apple by mediating the jaz1-trb1-myb9 complex. *Plant J.* 106, 1414–1430. doi: 10.1111/tpj.15245

An, X.-H., Tian, Y., Chen, K.-Q., Wang, X.-F., and Hao, Y.-J. (2012). The apple wd40 protein mdttg1 interacts with bhlh but not myb proteins to regulate anthocyanin accumulation. *J. Plant Physiol.* 169, 710–717. doi: 10.1016/j.jplph.2012.01.015

Andrews, S. (2010). *Fastqc: a quality control tool for high throughput sequence data.*

Bernardes, W. S., and Menossi, M. (2020). Plant 3' regulatory regions from mrna-encoding genes and their uses to modulate expression. *Front. Plant Sci.* 11. doi: 10.3389/fpls.2020.01252

Bianco, L., Cestaro, A., Sargent, D. J., Banchi, E., Derdak, S., Di Guardo, M., et al. (2014). Development and validation of a 20k single nucleotide polymorphism (snp)

## Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

## Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2024.1462545/full#supplementary-material

### SUPPLEMENTARY FIGURE 1

Heatmap-cluster analysis of differentially expressed genes (DEGs) between RF and WF. The standardized amount of expression is represented by different colors. GO-term and enrichment analysis of DEGs.

### SUPPLEMENTARY FIGURE 2

GO-enrichment analysis of differentially expressed genes (DEGs) between RF and WF GO terms were divided into three functional categories: molecular function, cellular component, and biological process. The Y-axis represents the number of DEGs enriched in each GO term.

### SUPPLEMENTARY FIGURE 3

Analysis of differentially methylated regions (DMRs) between RF and WF (A) Number of DMRs among the three contexts (CG, CHH, CHG). (B) Distribution of percentage of methylation for DMRs among the three contexts. (C) Position of DMRs associated with MdGSTU22 and MdPAL1.

whole genome genotyping array for apple (malus × domestica borkh). *PloS One* 9, e110377. doi: 10.1371/journal.pone.0110377

Bink, M. C. A. M., Jansen, J., Madduri, M., Voorrips, R. E., Durel, C.-E., Kouassi, A. B., et al. (2014). Bayesian qtl analyses using pedigreed families of an outcrossing species, with application to fruit firmness in apple. *Theor. Appl. Genet.* 127, 1073–1090. doi: 10.1007/s00122-014-2281-3

Bouillon, P., Fanciullino, A.-L., Belin, E., Bréard, D., Boisard, S., Bonnet, B., et al. (2024a). Image analysis and polyphenol profiling unveil red-flesh apple phenotype complexity. *Plant Methods* 20. doi: 10.1186/s13007-024-01196-1

Bouillon, P., Fanciullino, A.-L., Belin, E., Hanteville, S., Muranty, H., Bernard, F., et al. (2024b). Tracing the color: Qtl analysis reveals new insights into red-flesh pigmentation in apple (malus domestica). *Horticulture Res.* doi: 10.1093/hr/uhae171

Brueggemann, J., Weisshaar, B., and Sagasser, M. (2010). A wd40-repeat gene from malus × domestica is a functional homologue of arabidopsis thaliana transparent testa glabra1. *Plant Cell Rep.* 29, 285–294. doi: 10.1007/s00299-010-0821-0

Buhrman, K., Aravena-Calvo, J., Ross Zaulich, C., Hinz, K., and Laursen, T. (2022). Anthocyanic vacuolar inclusions: From biosynthesis to storage and possible applications. *Front. Chem.* 10. doi: 10.3389/fchem.2022.913324

Castillejo, C., Waurich, V., Wagner, H., Ramos, R., Oiza, N., Muñoz, P., et al. (2020). Allelic variation of myb10 is the major force controlling natural variation in skin and
flesh color in strawberry (fragaria spp.) fruit. *Plant Cell* 32, 3723-3749. doi: 10.1105/
tpc.20.00474

Celton, J., Gaillard, S., Bruneau, M., Pelletier, S., Aubourg, S., Martin-Magniette, M.,
et al. (2014). Widespread anti-sense transcription in apple is correlated with rna
production and indicates a large potential for transcriptional and/or post-
transcriptional control. *New Phytol.* 203, 287-299. doi: 10.1111/nph.12787

Chagné, D., Carlisle, C. M., Blond, C., Volz, R. K., Whitworth, C. J., Oraguzie, N. C.,
et al. (2007). Mapping a candidate gene (mdmyb10) for red flesh and foliage colour in
apple. *BMC Genomics* 8, 212. doi: 10.1186/1471-2164-8-212

Chagné, D., Lin-Wang, K., Espley, R. V., Volz, R. K., How, N. M., Rouse, S., et al. (2012). An ancient duplication of apple myb transcription factors is responsible for novel red fruit-flesh phenotypes. *Plant Physiol.* 161, 225-239. doi: 10.1104/
pp.112.206771

Chen, Z., Yu, L., Liu, W., Zhang, J., Wang, N., and Chen, X. (2021). Research progress of fruit color development in apple (malus domestica borkh.). *Plant Physiol. Biochem.* 162, 267-279. doi: 10.1016/j.plaphy.2021.02.033

Cheng, C.-Y., Krishnakumar, V., Chan, A. P., Thibaud-Nissen, F., Schobel, S., and
Town, C. D. (2017). enAraport11: a complete reannotation of the Arabidopsis thaliana
reference genome. *Plant J.* 89, 789-804. doi: 10.1111/tpj.13415

Daccord, N., Celton, J.-M., Linsmith, G., Becker, C., Choise, N., Schijlen, E.,
et al. (2017). High-quality *de novo* assembly of the apple genome and
methylene dynamics of early fruit development. *Nat. Genet.* 49, 1099-1106.
doi: 10.1038/ng.3886

Davies, K. M., Albert, N. W., and Schwinn, K. E. (2012). From landing lights to mimicry: the molecular regulation of flower colouration and mechanisms for pigmentation patterning. *Funct. Plant Biol.* 39, 619. doi: 10.1071/fp12195

Denoyes, B., Prohaska, A., Petit, J., and Rothan, C. (2023). Deciphering the genetic architecture of fruit color in strawberry. *J. Exp. Bot.* 74, 6306-6320. doi: 10.1093/jxb/
erad245

Ding, B., Patterson, E. L., Holalu, S. V., Li, J., Johnson, G. A., Stanley, L. E., et al. (2020). Two myb proteins in a self-organizing activator-inhibitor system produce spotted pigmentation patterns. *Curr. Biol.* 30, 802–814.e8. doi: 10.1016/
j.cub.2019.12.067

Eichenberger, M., Schwander, T., Hüppi, S., Kreuzer, J., Mittl, P. R. E., Peccati, F.,
et al. (2023). The catalytic role of glutathione transferases in heterologous anthocyanin
biosynthesis. *Nat. Catalysis* 6, 927–938. doi: 10.1038/s41929-023-01018-y

Espley, R. V., Bovy, A., Bava, C., Jaeger, S. R., Tomes, S., Norling, C., et al. (2012). Analysis of genetically modified red-fleshed apples reveals effects on growth and consumer attributes. *Plant Biotechnol. J.* 11, 408–419. doi: 10.1111/pbi.12017

Espley, R. V., Brendolise, C., Chagne, D., Kutty-Amma, S., Green, S., Volz, R., et al. (2009). Multiple repeats of a promoter segment causes transcription factor autoregulation in red apples. *Plant Cell* 21, 168–183. doi: 10.1105/tpc.108.059329

Espley, R. V., Hellens, R. P., Putterill, J., Stevenson, D. E., Kutty-Amma, S., and Allan, A. C. (2007). Red colouration in apple fruit is due to the activity of the myb transcription factor, mdmyb10. *Plant J.* 49, 414–427. doi: 10.1111/j.1365-313x.2006.02964.x

Espley, R. V., and Jaakola, L. (2023). The role of environmental stress in fruit
pigmentation. *Plant Cell Environ.* 46, 3663-3679. doi: 10.1111/pce.14684

Fairnie, A., Yeo, M., Gatti, S., Chan, E., Travaglia, V., Walker, J., et al. (2022). Eco-
evo-devo of petal pigmentation patterning. *Essays Biochem.* 66, 753–768. doi: 10.1042/
ebc20220051

Feng, G., Burleigh, J. G., Braun, E. L., Mei, W., and Barbazuk, W. B. (2017). Evolution of the 3r-myb gene family in plants. *Genome Biol. Evol.* 9, 1013–1029. doi: 10.1093/gbe/
evx056

Gonzalez, A., Zhao, M., Leavitt, J. M., and Lloyd, A. M. (2007). Regulation of the anthocyanin biosynthetic pathway by the ttg1/bhlh/myb transcriptional complex in arabidopsis seedlings. *Plant J.* 53, 814–827. doi: 10.1111/j.1365-313x.2007.03373.x

Grotewold, E. (2006). The genetics and biochemistry of floral pigments. *Annu. Rev.
Plant Biol.* 57, 761-780. doi: 10.1146/annurev.arplant.57.032905.105248

Han, M.-L., Yin, J., Zhao, Y.-H., Sun, X.-W., Meng, J.-X., Zhou, J., et al. (2020). How the color fades from malus halliana flowers: Transcriptome sequencing and dna methylation analysis. *Front. Plant Sci.* 11. doi: 10.3389/fpls.2020.576054

Honda, C., and Moriya, S. (2018). Anthocyanin biosynthesis in apple fruit.
*Horticulture* J. 87, 305-314. doi: 10.2503/hortj.okd-r01

Howard, N. P., Troggio, M., Durel, C.-E., Muranty, H., Denancé, C., Bianco, L., et al. (2021). Integration of infinium and axiom snp array data in the outcrossing species malus × domestica and causes for seemingly incompatible calls. *BMC Genomics* 22. doi: 10.1186/s12864-021-07565-7

Huang, L., Lin, B., Hao, P., Yi, K., Li, X., and Hua, S. (2024). Multi-omics analysis reveals that anthocyanin degradation and phytohormone changes regulate red color fading in rapeseed (brassica napus l.) petals. *Int. J. Mol. Sci.* 25, 2577. doi: 10.3390/
ijms25052577

Inc., P. T. (2015). *Collaborative data science.*

Jia, D., Wu, P., Shen, F., Li, W., Zheng, X., Wang, Y., et al. (2021). Genetic variation in the promoter of an r2r3myb transcription factor determines fruit malate content in apple (malus domestica borkh.). *Plant Physiol.* 186, 549–568. doi: 10.1093/plphys/
kiab098

Jiang, S.-H., Sun, Q.-G., Chen, M., Wang, N., Xu, H.-F., Fang, H.-C., et al. (2019).
Methylome and transcriptome analyses of apple fruit somatic mutations reveal the
difference of red phenotype. *BMC Genomics* 20. doi: 10.1186/s12864-019-5499-2

Jung, S., Lee, T., Cheng, C.-H., Buble, K., Zheng, P., Yu, J., et al. (2018). 15 years of
gdr: New data and functionality in the genome database for rosaceae. *Nucleic Acids Res.*
47, D1137–D1145. doi: 10.1093/nar/gky1000

Kim, D., Jeon, S. J., Yanders, S., Park, S., Kim, H. S., and Kim, S. (2022). Myb3 plays an important role in lignin and anthocyanin biosynthesis under salt stress condition in arabidopsis. *Plant Cell Rep.* 41, 1549–1560. doi: 10.1007/s00299-022-02878-7

Kokalj, D., Zlatić, E., Gićić, B., Kobav, M. B., and Vidrih, R. (2019). Postharvest
flavonol and anthocyanin accumulation in three apple cultivars in response to blue-
light-emitting diode light. *Scientia Hortic.* 257, 108711. doi: 10.1016/
j.scienta.2019.108711

Kondo, S., and Miura, T. (2010). Reaction-diffusion model as a framework for
understanding biological pattern formation. *Science* 329, 1616–1620. doi: 10.1126/
science.1179047

Kostick, S. A., and Luby, J. J. (2022). Apple fruit size qtls on chromosomes 8 and 16 characterized in 'honeycrisp'-derived germplasm. *Agronomy* 12, 1279. doi: 10.3390/
agronomy12061279

Krueger, F., and Andrews, S. R. (2011). Bismark: a flexible aligner and methylation caller for bisulfite-seq applications. *Bioinformatics* 27, 1571–1572. doi: 10.1093/
bioinformatics/btr167

Kumar, S., Deng, C. H., Molloy, C., Kirk, C., Plunkett, B., Lin-Wang, K., et al. (2022).
Extreme-phenotype gwas unravels a complex nexus between apple (It;igt;malus
domestica;igt); red-flesh colour and internal flesh browning. *Fruit Res.* 2, 1–14.
doi: 10.4813/frures-2022-0012

LaFontain, A. M., and Yuan, Y. (2021). Repressors of anthocyanin biosynthesis.
*New Phytol.* 231, 933-949. doi: 10.1111/nph.17397

Lallemand, T., Leduc, M., Desmazières, A., Aubourg, S., Rizzon, C., Landès, C., et al. (2023). Insights into the evolution of ohnologous sequences and their epigenetic marks post-wgd in malus domestica. *Genome Biol. Evol.* 15. doi: 10.1093/gbe/evd178

Lev-Yadun, S. (2016). "Ant Mimicry," in *Defensive (anti-herbivory) Coloration in Land Plants* (Springer, Cham), p. 299-304. doi: 10.1007/978-3-319-42096-7_58

Li, H., Duan, S., Sun, W., Wang, S., Zhang, J., Song, T., et al. (2022). Identification, through transcriptome analysis, of transcription factors that regulate anthocyanin biosynthesis in different parts of red-fleshed apple 'may' fruit. *Hortic. Plant J.* 8, 11–21. doi: 10.1016/j.hpj.2021.07.001

Li, M., Coneva, V., Robbins, K. R., Clark, D., Chitwood, D., and Frank, M. (2021).
Quantitative dissection of color patterning in the foliar ornamental coleus. *Plant
Physiol.* 187, 1310-1324. doi: 10.1093/plphys/kiab393

Li, Z., Liu, W., Chen, Q., Zhang, S., Mei, Z., Yu, L., et al. (2023). Mdm-mir858 targets
mdmyb9 and mdmybpa1 to participate anthocyanin biosynthesis in red-fleshed apple.
*Plant J.* 113, 1295-1309. doi: 10.1111/tpj.16111

Liu, H., Liu, Z., Wu, Y., Zheng, L., and Zhang, G. (2021). Regulatory mechanisms of
anthocyanin biosynthesis in apple and pear. *Int. J. Mol. Sci.* 22, 8441. doi: 10.3390/
ijms22168441

Li, W., Wang, Y., Yu, L., Jiang, H., Guo, Z., Xu, H., et al. (2019). Mdwrkyl1
participates in anthocyanin accumulation in red-fleshed apples by affecting myb
transcription factors and the photoresponse factor mdhy5. *J. Agric. Food Chem.* 67,
8783-8793. doi: 10.1021/acs.jafc.9b02920

Lloyd, A., Brockman, A., Aguirre, L., Campbell, A., Bean, A., Cantero, A., et al. (2017). Advances in the myb-bhlh-wd repeat (mbw) pigment regulatory model: Addition of a wrky factor and cooption of an anthocyanin myb for betalain regulation. *Plant Cell Physiol.* 58, 1431–1441. doi: 10.1093/pcp/pcx075

Love, M. I., Huber, W., and Anders, S. (2014). Moderated estimation of fold change and dispersion for rna-seq data with deseq2. *Genome Biol.* 15. doi: 10.1186/s13059-014-0550-8

Lu, Z., Cao, H., Pan, L., Niu, L., Wei, B., Cui, G., et al. (2021). Two loss-of-function alleles of the glutathione s-transferase (gst) gene cause anthocyanin deficiency in flower and fruit skin of peach (prunus persica). *Plant J.* 107, 1320–1331. doi: 10.1111/tpj.15312

Meng, J.-X., Wei, J., Chi, R.-F., Qiao, Y.-H., Zhou, J., Wang, Y.-L., et al. (2022).
Mrmby44-like negatively regulates anthocyanin biosynthesis and causes spring leaf
color of malus 'radiant' to fade from red to green. *Front. Plant Sci.* 13. doi: 10.3389/
fpbs.2022.822340

Moyers, B. T., Owens, G. L., Baute, G. J., and Rieseberg, L. H. (2017). The genetic
architecture of uv floral patterning in sunflower. *Ann. Bot.* 120, 39–50. doi: 10.1093/
aob/mcx038

Pelletier, S. (2022). *Anadiff: A tool for differential analysis of microarrays and rnaseq.*
doi: 10.5281/ZENODO.6477918

Piepho, H. P., Möhring, J., Melchinger, A. E., and Büchse, A. (2007). Blup for phenotypic selection in plant breeding and variety testing. *Euphytica* 161, 209–228. doi: 10.1007/s10681-007-9449-8

Porebski, S., Bailey, L. G., and Baum, B. R. (1997). Modification of a ctab dna extraction protocol for plants containing high polysaccharide and polyphenol components. *Plant Mol. Biol. Rep.* 15, 8–15. doi: 10.1007/bf02772108

Powell, A. A., Kostick, S. A., Howard, N. P., and Luby, J. J. (2022). Elucidation and characterization of qtls for russet formation on apple fruit in 'honeycrisp'-derived breeding germplasm. *Tree Genet. Genomes* 19. doi: 10.1007/s11295-022-01582-7
Quevillon, E., Silventoinen, V., Pillai, S., Harte, N., Mulder, N., Apweiler, R., et al. (2005). InterProScan: protein domains identifier. *Nucleic Acids Res.* 33, W116–W120. doi: 10.1093/nar/gki442

Quinlan, A. R., and Hall, I. M. (2010). Bedtools: a flexible suite of utilities for comparing genomic features. *Bioinformatics* 26, 841–842. doi: 10.1093/bioinformatics/btq033

Rong, H., Han, X., Xin, Y., Ni, Z., Zhang, W., and Xu, L. (2023). Small rna and degradome sequencing reveal roles of mirnas in the petal color fading of malus crabapple. *Int. J. Mol. Sci.* 24, 11384. doi: 10.3390/ijms241411384

Shen, F., Hu, C., Huang, X., Wu, R., Luo, S., Xu, C., et al. (2022). Characterization of the genetic and regulatory networks associated with sugar and acid metabolism in apples via an integrated strategy. *Front. Plant Sci.* 13. doi: 10.3389/fpls.2022.1066592

Shi, L., Li, X., Fu, Y., and Li, C. (2023). Environmental stimuli and phytohormones in anthocyanin biosynthesis: A comprehensive review. *Int. J. Mol. Sci.* 24, 16415. doi: 10.3390/ijms242216415

Slatnar, A., Licznar-Malanczuk, M., Mikulic-Petkovsek, M., Stampar, F., and Veberic, R. (2014). Long-term experiment with orchard floor management systems: Influence on apple yield and chemical composition. *J. Agric. Food Chem.* 62, 4095–4103. doi: 10.1021/jf500846r

Tanaka, Y., Sasaki, N., and Ohmiya, A. (2008). Biosynthesis of plant pigments: anthocyanins, betalains and carotenoids. *Plant J.* 54, 733–749. doi: 10.1111/j.1365-313x.2008.03447.x

Tzelias, A., Lin-Wang, K., Stevenson, D. E., Cooney, J. M., Hellens, R. P., Allan, A. C., et al. (2011). Apple skin patterning is associated with differential expression of myb10. *BMC Plant Biol.* 11, 93. doi: 10.1186/1471-2229-11-93

Todesco, M., Bercovich, N., Kim, A., Imerovski, I., Owens, G. L., Ruiz, D., et al. (2022). Genetic basis and dual adaptive role of floral pigmentation in sunflowers. *eLife* 11. doi: 10.7554/elife.72072

Van der Walt, S., Schönberger, J. L., Nunez-Iglesias, J., Boulogne, F., Warner, J. D., Yager, N., et al. (2014). scikit-image: image processing in python. *PeerJ* 2, e453. doi: 10.7717/peerj.453

Vanderzande, S., Howard, N. P., Cai, L., Da Silva Linge, C., Antanaviciute, L., Bink, M. C. A. M., et al. (2019). High-quality, genome-wide snp genotypic data for pedigree germplasm of the diploid outbreeding species apple, peach, and sweet cherry through a common workflow. *PloS One* 14, e0210928. doi: 10.1371/journal.pone.0210928

van Nocker, S., Berry, G., Najdowski, J., Michelutti, R., Luffman, M., Forsline, P., et al. (2011). Genetic diversity of red-fleshed apples (malus). *Euphytica* 185, 281–293. doi: 10.1007/s10681-011-0579-7

Verma, S., Evans, K., Guan, Y., Luby, J. J., Rosyara, U. R., Howard, N. P., et al. (2019). Two large-effect qtls, ma and ma3, determine genetic potential for acidity in apple fruit: breeding insights from a multi-family study. *Tree Genet. Genomes* 15. doi: 10.1007/s11295-019-1324-y

Verweij, W., Spelt, C. E., Bliek, M., de Vries, M., Wit, N., Faraco, M., et al. (2016). Functionally similar wrky proteins regulate vacuolar acidification in petunia and hair development in arabidopsis. *Plant Cell* 28, 786–803. doi: 10.1105/tpc.15.00608

Volz, R. K., McGhie, T. K., and Kumar, S. (2014). Variation and genetic parameters of fruit colour and polyphenol composition in an apple seedling population segregating for red leaf. *Tree Genet. Genomes* 10, 953–964. doi: 10.1007/s11295-014-0734-0

Voorrips, R. E., Bink, M. C. A. M., Kruisselbrink, J. W., Koehorst-van Putten, H. J. J., and van de Weg, W. E. (2016). Pediplaphotyper: software for consistent assignment of marker haplotypes in pedigrees. *Mol. Breed.* 36. doi: 10.1007/s11032-016-0539-y

Wang, N., and Chen, X. (2021). “Genetics and Genomics of Fruit Color Development in Apple,” in *The Apple Genome. Compendium of Plant Genomes*. Ed. Korban, S. S. (Springer, Cham). doi: 10.1007/978-3-030-74682-7_13

Wang, N., Liu, W., Zhang, T., Jiang, S., Xu, H., Wang, Y., et al. (2018). Transcriptomic analysis of red-fleshed apples reveals the novel role of mdwrky11 in flavonoid and anthocyanin biosynthesis. *J. Agric. Food Chem.* 66, 7076–7086. doi: 10.1021/acs.jafc.8b01273

Wang, W., Moss, S. M. A., Zeng, L., Espley, R. V., Wang, T., Lin-Wang, K., et al. (2022). The red flesh of kiwifruit is differentially controlled by specific activation-repression systems. *New Phytol.* 235, 630–645. doi: 10.1111/nph.18122

Wang, N., Zheng, Y., Duan, N., Zhang, Z., Ji, X., Jiang, S., et al. (2015). Comparative transcriptomes analysis of red- and white-fleshed apples in an f1 population of malus sieversii f. niedzwetzkyana crossed with m. domestica 'fuji'. *PloS One* 10, e0133468. doi: 10.1371/journal.pone.0133468

Wei, L., Mao, W., Jia, M., Xing, S., Ali, U., Zhao, Y., et al. (2018). Famyb44.2, a transcriptional repressor, negatively regulates sucrose accumulation in strawberry receptacles through interplay with famyb10. *J. Exp. Bot.* 69, 4805–4820. doi: 10.1093/jxb/ery249

Xiang, Y., Yu Huang, X., Zhao, Y.-W., Wang, C.-K., Sun, Q., and Hu, D.-G. (2024). Role of an atp binding cassette (abc) transporter mdabc17 in the anthocyanin accumulation of apple. *Scientia Hortic.* 323, 112502. doi: 10.1016/j.scienta.2023.112502

Xu, S. (2022). "Estimation of Heritability," in *Quantitative Genetics* (Springer, Cham), p. 147–175. doi: 10.1007/978-3-030-83940-6_10

Xu, J., Xiong, L., Yao, J.-L., Zhao, P., Jiang, S., Sun, X., et al. (2024). Hypermethylation in the promoter regions of flavonoid pathway genes is associated with skin color fading during 'daihong' apple fruit development. *Horticulture Res.* 11. doi: 10.1093/hr/uhae031

Xu, H., Zou, Q., Yang, G., Jiang, S., Fang, H., Wang, Y., et al. (2020). Mdmyb6 regulates anthocyanin formation in apple both through direct inhibition of the biosynthesis pathway and through substrate removal. *Horticulture Res.* 7. doi: 10.1038/s41438-020-0294-4

Yan, H., Pei, X., Zhang, H., Li, X., Zhang, X., Zhao, M., et al. (2021). Myb-mediated regulation of anthocyanin biosynthesis. *Int. J. Mol. Sci.* 22, 3103. doi: 10.3390/ijms22063103

Yuan, Y. (2018). Monkeyflowers (mimulus): new model for plant developmental genetics and evo-devo. *New Phytol.* 222, 694–700. doi: 10.1111/nph.15560

Zhang, Y., Butelli, E., and Martin, C. (2014). Engineering anthocyanin biosynthesis in plants. *Curr. Opin. Plant Biol.* 19, 81–90. doi: 10.1016/j.pbi.2014.05.011

Zhang, S., Chen, Y., Zhao, L., Li, C., Yu, J., Li, T., et al. (2020). A novel nac transcription factor, mdnac42, regulates anthocyanin accumulation in red-fleshed apple by interacting with mdmyb10. *Tree Physiol.* 40, 413–423. doi: 10.1093/ treephys/tpaa004

Zhang, H., Lang, Z., and Zhu, J.-K. (2018). Dynamics and function of dna methylation in plants. *Nat. Rev. Mol. Cell Biol.* 19, 489–506. doi: 10.1038/s41580- 018-0016-z

Zhang, W., Liao, L., Xu, J., Han, Y., and Li, L. (2021). Genome-wide identification, characterization and expression analysis of mate family genes in apple (malus × domestica borkh). *BMC Genomics* 22. doi: 10.1186/s12864-021-07943-1

Zhang, X., Zhang, K., Guo, Y., Lv, X., Wu, M., Deng, H., et al. (2023). Methylation of acgt1 is associated with anthocyanin accumulation in the kiwifruit outer pericarp. *J. Agric. Food Chem.* 71, 18865–18876. doi: 10.1021/acs.jafc.3c03029

Zhao, Y.-W., Wang, C.-K., Huang, X.-Y., and Hu, D.-G. (2021b). Genome-wide analysis of the glutathione s-transferase (gst) genes and functional identification of mdgstu12 reveals the involvement in the regulation of anthocyanin accumulation in apple. *Genes* 12, 1733. doi: 10.3390/genes12111733

Zhao, G., Xiang, F., Zhang, S., Song, J., Li, X., Song, L., et al. (2021a). Pblac4-like, activated by pbmyb26, related to the degradation of anthocyanin during color fading in pear. *BMC Plant Biol.* 21. doi: 10.1186/s12870-021-03220-1

Zheng, X., Om, K., Stanton, K. A., Thomas, D., Cheng, P. A., Eggert, A., et al. (2021). The regulatory network for petal anthocyanin pigmentation is shaped by the myb5a/negan transcription factor in mimulus. *Genetics* 217. doi: 10.1093/genetics/iyaa036

Zhou, H., Lin-Wang, K., Wang, F., Espley, R. V., Ren, F., Zhao, J., et al. (2018). Activator-type r2r3-myb genes induce a repressor-type r2r3-myb gene to balance anthocyanin and proanthocyanidin accumulation. *New Phytol.* 221, 1919–1934. doi: 10.1111/nph.15486

Zhou, M., Zhang, K., Sun, Z., Yan, M., Chen, C., Zhang, X., et al. (2017). Lnk1 and Lnk2 corepressors interact with the myb3 transcription factor in phenylpropanoid biosynthesis. *Plant Physiol.* 174, 1348–1358. doi: 10.1104/pp.17.00160

Zingaretti, L. M., Monfort, A., and Pérez-Enciso, M. (2021). Automatic fruit morphology phenotype and genetic analysis: An application in the octoploid strawberry. *Plant Phenomics* 2021. doi: 10.34133/2021/9812910