Prediction of fruit texture with training population optimization for efficient genomic selection in apple

Texture plays a major role in the determination of fruit quality in apple. Due to its physiological and economic relevance, this trait has been largely investigated, leading to the fixation of the major gene PG1 controlling firmness in elite cultivars. To further improve fruit texture, the targeting of an undisclosed reservoir of loci with minor effects is compelling. In this work, we aimed to unlock this potential with a genomic selection approach by predicting fruit acoustic and mechanical features as obtained with a TA.XTplus texture analyzer in 537 individuals genotyped with 8,294 SNP markers. The best prediction accuracies following cross-validations within the training set (TRS) of 259 individuals were obtained for the acoustic linear distance (0.64). Prediction accuracy was further improved through the optimization of TRS size and composition according to the test set. With this strategy, a maximal accuracy of 0.81 was obtained when predicting the synthetic trait PC1 in the family ‘Gala × Pink Lady’. We discuss the impact of genetic relatedness and clustering on trait variability and predictability. Moreover, we demonstrated the need for a comprehensive dissection of the complex texture phenotype and the potentiality of using genomic selection to improve fruit quality in apple. Highlight A genomic selection study, together with the optimization of the training set, demonstrated the possibility to accurately predict texture sub-traits valuable for the amelioration of fruit quality in apple.

Fruits, during maturation and ripening, undergo a complex series of genetically 65 programmed events contributing to their attractiveness and suitability for human 66 consumption. Amongst the various physiological and physical changes, fruit texture is 67 certainly the most important and investigated traits, especially in apple. A favorable texture is 68 in fact highly appreciated by consumers, enabling, moreover, a long-term storage. 69 Texture can nowadays be dissected into two groups of sub-traits, mechanical and 70 acoustic, contributing to distinguish between firm (based on mechanical sub-traits) and crispy 71 (based on acoustic sub-traits) types of apples. These texture parameters have been already 72 described and validated in apple (Costa et al., 2011(Costa et al., , 2012, and were implemented in QTL-73 mapping studies carried out with bi-parental populations (Longhi et al., 2012) as well as more 74 structured approaches, such as Pedigreed Based Analysis (PBA) and Genome-Wide 75 Association Studies (GWAS, Kumar et al., 2013;Migicovsky et al., 2016;Amyotte et al., 76 2017;Di Guardo et al., 2017;McClure et al., 2019). These works elucidated the complex 77 genetic control of the fruit texture in apple, identifying a large number of QTLs distributed 78 over the apple genome, with the most relevant regions located on chromosome 3, 10 and 16. 79 This genetic complexity is moreover reflected in the regulation of the cell-wall and middle 80 grown according to conventional horticultural management for plant training, pruning and 149 pest-disease control. 150 Fruits were harvested from each plant at the time of the physiological ripening stage, 151 established according to standard horticultural fruit quality parameters, such as the change in 152 color of the skin, seeds and flesh, fruit firmness value and the iodine coloration index 153 indicating the internal starch degradation. After harvest, fruits were stored for two months at 154 2°C with 95% of relative humidity. 155 156 Texture phenotyping 157 The texture performance of the apple fruit was phenotypically dissected into 158 mechanical and acoustic sub-traits with the use of a texture analyzer TA.XTplus (Stable 159 MicroSystems Ltd., Godalming, UK) equipped with an acoustic envelop device AED (Stable 160 MicroSystems Ltd., Godalming, UK), as described in Costa et al. (2011). For each genotype 161 included in the population, a homogeneous set of five apples was collected. Four identical 162 discs were isolated per fruit, avoiding seeds, seed cavity tissues or skin, for a total of 20 163 measurements per genotype (5 biological replicates and 4 technological replicates). Each 164 texture profile was then digitally elaborated identifying 12 texture measurements (i. e. 'sub-165 traits'), four related to the acoustic performance and eight to the mechanical force-166 displacement. In brief, the mechanical sub-traits were coded as: initial, final, maximum and 167 mean force (related to the different force values associated to the different parts of the force-168 displacement profile), area, force linear distance (derived length of the profile), Young's 169 module (also known as elasticity module) and number of force peaks. The acoustic sub-traits 170 were maximum and mean acoustic pressure, acoustic linear distance and number of acoustic 171 peaks. A more exhaustive and complete description of the texture sub-traits is reported in 172 Costa et al. 2011. 173 174 SNP genotyping 175 The DNA employed for the genotyping of each individual considered in this survey 176 was isolated from young leaves collected at the beginning of the vegetative phase with the 177 Qiagen DNeasy Plant Kit and further quantified with a Nanodrop ND-8000 178 (ThermoScientific, USA). SNP markers were genotyped through the HiScan (Illumina, USA) 179 and the apple 20K SNP chip Infinium array (Illumina, USA) assembled within the framework 180 of the European project FruitBreedomics (Bianco et al., 2014). The SNP pattern was initially 181 analyzed with the software GenomeStudio and further re-edited with ASSiST (Di Guardo et 182 al., 2015). SNPs with minor allele frequencies lower than 0.05 and call rate below 0.2 were 183 filtered out with the package 'snpStats' (Clayton, 2019). The final set of markers successfully 184 recovered in the population consisted in 8,294 biallelic SNPs. 185

Analysis of the fruit texture sub-traits 186
We used a mixed linear model to get the best linear unbiased predictors (BLUPs) of 187 each individual's genotypic value. For each apple measured, we first calculated the mean over 188 the four technical replicates to retain only the biological replication level in the model. Each 189 of the twelve mechanical or acoustic sub-traits, considered as 'Y', was explained by the 190 genotype as random effect, the trial (location by year) as fixed effect and the random effect of 191 the error as:  (Lê et al., 2008). Only values from the collection were used to create the principal 199 The realized additive relationship was calculated with the 'A.mat' function of the 206 'rrBLUP' package (Endelman, 2011) and depicted in a heatmap plot obtained with the R-207 function 'heatmap.2' (package 'gplots', Warnes et al., 2016). Genetic clustering was further 208 assessed in the collection with a discriminant analysis of principal components (DAPC, 209 Jombart et al., 2010), carried out with the R-package 'adegenet' (Jombart, 2008) using the 210 entire set of 8,294 markers. In the first step, six significant clusters were retained with the 211 function 'find.clusters' using 300 principal components and selecting the number of clusters 212 with the highest likelihood (based on the Bayesian information criterion value-BIC, Fig. S1). 213 Out of these variables, 150 were retained and employed in the clustering computed with the 214 'dapc' function, which created five principal components that maximized the inter-cluster 215 distance while minimizing the inter-individual distance within each cluster. The assignment of 216 offsprings to clusters was obtained with the function 'predict_dapc'. Pairwise Fst values 217 between clusters were then computed with the entire SNP set with the function 218 'pairwise.WCfst' from R-package 'hierfstat' (Yang, 1998, Goudet 2005. 219 220

Prediction models 221
Genomic predictions were computed through two models implemented in the rrBLUP 222 framework, as reported in Endelman et al. 2011 (and 'rrBLUP R'- A 5-fold cross-validation was applied within the collection with both model (A) and 234 (B) respectively and repeated 100 times. For predicting each family (considered then as TS), 235 three different TRS composition rules, named as "scenarios", were tested using the two 236 models without a priori genetic information on individuals. In scenario 1, each family was 237 predicted using the collection only. In scenario 2, 30% of individuals of the predicted family 238 were instead added to the collection in the TRS while the remaining 70% formed the TS. In 239 scenario 3, a single half-sib family (e.g. 'GaPL' is half-sib with 'FjPL' and 'GaPi') was added 240 to the collection to form the TRS, leading to two to four TRS possibilities (and accuracy 241 values). To illustrate the scenarios taking 'GaPi' as an example, scenario 1 corresponded to 242 To this end, a relatedness-driven and a principal component-driven approaches were adopted. 250 The relatedness-driven approach was tested in three different manners: (i) by starting with the 251 10 most-related individuals and adding single individuals with decreasing mean relationship 252 to the family, or (ii) with decreasing maximum relationship to the family (N=10 to N=259), or 253 (iii) by starting with a TRS composed of the most related cluster and adding less and less 254 related clusters successively (final TRS size N=259). In the principal component-driven 255 approach, TRS individuals were selected with increasing TRS size using a protocol by 256 Akdemir (R-package 'STPGA', 2019). The optimal TRS with increasing size from 10 257 individuals to 259 with increments of 20 individuals was chosen based on the five principal 258 components obtained with DAPC analysis and using the 'CDmean' design criteria and the 259 function 'GenAlgForSubsetSelection'. Here, individuals were chosen independently for each 260 TRS size, meaning that we did not proceed to a gradual enrichment of the TRS. 261 All accuracy values were based on Pearson correlation calculated between observed 262 values (i.e. BLUPs of genotypic values) and predicted values of the TS individuals. When 263 standard deviations were not available, we calculated an approximate 95% confidence interval 264 of the correlation coefficient with a Fisher's Z-transformation ('cor.test' function in base R). 265 Calculations were performed in R (R Core Team, 2014) and graphs were created with the R-266 package 'ggplot2' (Wickham, 2016). 267

269
Fruit texture phenotypic dissection 270 The fruit texture phenotypic data used in this survey were represented by the analysis 271 of multi-trait features accurately dissected into 4 acoustic and 8 mechanical sub-traits (Table  272 2, Table S1). A mixed linear model was used to obtain BLUPs of genotypic values used in the 273 further analyses. The texture sub-traits showed an overall high heritability, spanning from 274 0.90-0.96 for the entire population (collection and families) to 0.88-0.94 for the apple 275 accessions included in the collection (Table 2). In order to visualize the diversity and 276 inheritance of fruit texture profiles, a principal component analysis (PCA) was performed 277 using the twelve textural sub-traits measured in the collection, while individuals from families 278 were considered as supplementary individuals (see also Di Guardo et al. 2017, Fig. 1). In this 279 analysis, the first PC axis (PC1), explaining 80.5% of phenotypic variability, comprehensively 280 summarizing the general variability of the twelve phenotypic variables. The second axis 281 (PC2), instead, mainly differentiated the acoustic from mechanical sub-traits, explaining a 282 smaller, yet substantial, portion of the phenotypic variability (12.7%, Fig. 1A). 283 In the distinction between the two types of texture sub-traits (mechanical and acoustic) 284 by PC2, it is worth noting that one mechanical variable (FNP) was oriented together with the 285 acoustic group. FNP was in fact more correlated with acoustic sub-traits (mean correlation 286 0.77) than with the rest of the mechanical ones (mean correlation 0.69, Fig. 1A). Individuals 287 of the population were present in the four quadrants of the PCA 2D-plot, identifying different 288 types of texture: mealy (negative PC1), predominantly firm (positive PC1 and negative PC2) 289 and predominantly crispy (positive PC1 and positive PC2, Fig. 1B). With this regard, the 290 distribution of texture profiles indicated that the collection is mainly composed of individuals 291 with low to moderate crispiness and firmness at the exception of few outliers. It is also 292 important to note that variation on the PC2 axis is much lower for accessions having a 293 negative PC1 value, illustrating that mealy apples cannot be crispy (Fig. 1B). 294 The six parental cultivars, known to have different texture profiles after two months of 295 storage, were, as expected, plotted over the different quadrants of the PCA 2D-plot (Fig. 1C). 296 'Delearly' and 'Golden Delicious' were plotted in the area corresponding to the mealy type of 297 apple, while 'Royal Gala' was instead grouped with moderately firm apples. 'Fuji', 'Pink 298 Lady' and 'Pinova' were instead positioned in the positive quadrant for both PC1 and PC2,299 corresponding to the crispy type of apple. The populations originated by the controlled cross 300 of these varieties were also distributed over the PCA plot with specific orientations (Fig. 1B-301 C). In particular, 'FjPL' offsprings were mostly projected towards the 'firm quadrant', while 302 'GDFj' was more oriented in the 'crispy quadrant' (Fig. 1B). Moreover, the segregation of the 303 families was very variable with regard to their corresponding parental profiles (Fig. 1C). 304 While 'GDFj' was the only family showing a classic type of segregation (intermediate 305 between the parents), the distributions of the other families were more similar to one of the 306 two parents ('FjDe' and 'GaPi'), with a varying number of offsprings being of transgressive 307 type ('FjDe' ,'GaPL', 'FjPi' and 'FjPL'). In particular, while 'Fuji' and 'Pink Lady' showed a 308 very similar texture profile on PC1 (2.99 and 3.14 respectively), major differences were 309 observed on the PC2 (1.6 and 0.51 respectively, Fig. 1C, Table S1). Variation in the texture 310 performance of 'FjPL' offsprings was also observed on the PC2 axis, although with a much 311 broader variation with regards to 'Fuji' and 'Pink Lady'. Accordingly, apples of this family 312 were overall firm to very firm while having a very low to very high crispiness (Fig. 1C, Table  313 S1, Fig. S2). 314 315 Additive relationship and genetic clustering in the population 316 The accuracy of genomic prediction is highly correlated to the level of relatedness 317 between the training and the test sets (TRS and TS). To identify the overall patterns of 318 relatedness between families and the collection, a clustering analysis of all the individuals 319 based on their pairwise additive relationship was performed (Fig. 2). The parental cultivar 320 'Royal Gala' was found to be the most related to the rest of the collection (mean additive 321 relatedness -6.32E-4), while 'Fuji' was the most distantly related (mean additive relatedness -322 0.102, Table S2). Accordingly, 'Royal Gala'-related families were more closely related to the 323 collection respect to the four 'Fuji'-related families, plotted together on the top-right panel of 324 the heatmap (Fig. 2). Mean additive relationship values for each family reflected the patterns 325 observed on the heatmap, namely higher values for 'GaPi' and 'GaPL' (-0.021 to -0.020) and 326 lower for 'Fuji'-related families (-0.056 to -0.078, Table 1, Table S2). 327 To investigate the genetic structure of the collection and its impact on the prediction 328 accuracy, a discriminant analysis of principal component with the entire SNP set (8,294  329 SNPs) was performed. Through the BIC criteria, six clusters, described with five principal 330 components, were defined as the most probable (see Methods, Fig. S1). All parental cultivars 331 were assigned to cluster 5, except 'Fuji' that was grouped in cluster 2 (Fig. 3A, Table S3). Of 332 these clusters, cluster 5 resulted to be the largest (N=66), while the smallest was cluster 6 333 (n=25, Table 1, Table S3). The cluster assignment in families was predicted using the 334 principal components derived by the DAPC analysis carried out on the collection. Most of the 335 individuals were assigned to the parental clusters 2 and 5, while 8 individuals of 'FjDe' and 336 one of 'FjPi' were assigned to cluster 1 (Table 1,  (models A and B, respectively). In this context, PC1 and PC2 were also considered as traits, 355 leading in the end to 14 predicted traits (Fig. S4). Instead of improving predictions, the 356 inclusion of the clustering effect degraded accuracies for all traits, with a maximum accuracy 357 decrease of 0.02 for the mean force (FMean). The highest mean prediction was obtained for 358 the acoustic linear distance (ALD, Fig. S4) whereas the number of force 359 peaks yielded the second highest accuracy (FNP, Fig. S4, Table S5). 360 Moreover, while FNP yielded a relatively high accuracy as inferred from heritability (0.93, 361  In practice, families can be predicted with any available related genetic material that 367 has been genotyped and phenotyped. For this reason, three different scenarios of training 368 population design were tested, including or not individuals from the predicted family or from 369 a half-sib family (see Methods, "Prediction models"). The predictions in each of these 370 scenarios were calculated with the two prediction models (A and B, respectively depicted in 371 Fig. 4, Fig. S5). Without clustering, overall three families ('FjPi', 'GaPi' and 'GaPL') could 372 be predicted with moderate to high accuracies (accuracies ranging from 0.08 for PC2 in 373 'GaPi' to 0.73 for PC1 in 'GaPL', respectively), with PC1 being the best predicted trait 374 among these families (mean for scenario 1, model A: 0.50, Fig. 4) Table S6, Fig 4, Fig. S5). Thus, the implementation of 382 clustering did not clearly improve the predictions of families. 383 It is also important to underline that the addition of related individuals to the collection 384 did not systematically improve the predictions. For instance, in 'GaPL' the prediction was 385 more accurate with scenario 1 with regards to scenario 2 and 3 (mean prediction accuracies of 386 0.60, 0.56 and 0.53 respectively for scenario 1, 2, 3, respectively, model A). Scenario 2 387 particularly improved the accuracies in 'FjPi' (mean accuracies of 0.32, model A) as it better 388 predicted 12 out of 14 traits. Scenario 3 instead was the lowest performing, although it 389 increased the prediction accuracy of 7 traits (8 with clustering) in 'GaPi' (mean accuracy of 390 0.38, all values across trait in model A, Fig 4, Table S6). 391 392

Genomic prediction of families with training population optimization 393
To test the hypothesis that retaining only the most related individuals or clusters in the 394 TRS might allow to maximize prediction accuracies, we compared the predictive abilities 395 obtained for each family and trait using training sets with different sizes. This process started 396 with a small TRS having the highest relatedness to which individuals were added in the order 397 of decreasing relatedness to reach the size of the entire collection using three different 398 enrichment procedures (see Methods). TRS optimization was also carried out with a more 399 sophisticated approach based on the optimization algorithm presented by Akdemir et al.  Table S7. Regarding the four selected traits, the best accuracy for each of the 405 6 ൈ 4 family-trait combinations was in most cases obtained with the addition of single 406 individuals based on their relationship to the family (in 10 cases using the maximum 407 relationship and in 10 cases using the mean relationship, Fig. 5A and B, Table 3, Table S7). 408 The mean optimal population size was 92 individuals with a minimum size of 10 and a 409 maximum size of 202 individuals (Table 3, Table S7), meaning that the entire collection was 410 never considered as the optimal TRS for predicting texture. The maximal accuracies observed 411 ranged from 0.01 to 0.81, which corresponded to a mean increase in accuracy of 0.17 when 412 compared to predictions of families with the entire collection (minimum increase: 0.02; 413 maximum increase: 0.40 -compared to scenario 1, model A). The highest accuracy was 0.81, 414 and was obtained for the "multi-trait" PC1 in 'GaPL' family with only 129 individuals, i.e. 415 nearly half of the collection size. The distribution of accuracies with increasing TRS size in 416 each family for the four focal traits was also investigated (Fig. 5). Overall, traits tended to 417 follow the same trend within a family. In families 'GaPL' and 'GaPi', which had the highest 418 relatedness to the collection among all families (Table 1), the accuracy was moderate to high 419 from as few as 100 individuals for ALD, FNP and PC1, and remained relatively stable while 420 increasing TRS size (Fig. 5A-D). 'FjPi' was the only family for which increasing TRS up to 421 200 individuals resulted in a clear accuracy improvement, with any of the approaches 422 implemented here (Fig. 5A-D). In families with overall low accuracies, such as 'FjDe', 'FjPL' 423 and 'GDFj', the highest accuracy was in most cases obtained with 10 to 70 individuals, and 424 declined or remained stable with larger TRS size ( Fig. 5A-D). In 'GDFj', for instance, 425 accuracies above 0.2 were found only with a TRS of 10 to 66 individuals ( Fig. 5A ), an improved accuracy of 0.32 was observed with as few as 15 individuals 428 (based on maximum relationship, Fig. 5B). 429

Discussion 430
In this work we assessed the feasibility of genomic selection (GS) for apple texture by 431 performing an in-depth analysis of this complex phenotype together with the genetic 432 correlates influencing its genomic predictions. The results presented here on genomic 433 prediction for apple texture evidenced a large potential for GS for this trait, providing 434 important key elements and tools to set-up a prediction experiment given the available genetic 435 information in any apple population. 436

437
Family-dependent fruit texture profiles and fruit texture prediction 438 The texture dissected "sub-traits" were highly heritable, although variability within 439 families was very contrasted, showing, in specific cases, a transgressive segregation, such as 440 'FjPL'. Although the traits were predictable with moderate to high accuracy within the 441 collection (accuracies between 0.41 and 0.64), this was not easily achievable in all biparental 442 families. Without TRS optimization, texture could be accurately predicted for 'GaPL' (mean 443 accuracy of 0.57), while 'GaPi' and in 'FjPi' showed a moderate prediction accuracy (mean 444 accuracy of 0.30). In contrast, near-zero or negative accuracies were instead obtained for 445 'FjDe', 'FjPL' and 'GDFj' across traits (mean accuracy of -0.05). Surprisingly, large negative 446 accuracy values were repeatedly obtained in 'FjDe' and 'GDFj', which could be potentially 447 explained by the strong epistatic effect possibly present in these families (Lehner, 2011) or by 448 a systematic bias due to the calculation of the Pearson correlation coefficient (Zhou et al., 449 2016), indicating that fruit texture cannot be predicted in these families using the entire 450 collection as TRS. In contrast, previous works on firmness and crispiness yielded mostly low 451 accuracies when predicting unobserved genotypes in a set of families or in a collection 452 (between 0.15 and 0.35, Kumar et al. 2015, McClure et al. 2018. A much higher accuracy of 453 0.83 was found for firmness by Kumar et al. (2012), which can be mainly explained by their 454 crossing design and validation procedure. In the present study, the analysis of PCA allowed to 455 better understand the relation between firmness and crispiness, both positively correlated and 456 summarized by PC1 and PC2, with PC2 specifically dissecting the difference between these 457 two texture sub-traits. When used as synthetic trait in the computation, PC1 was among the 458 best predictable traits (accuracy of 0.59 in collection and highest accuracy among traits and 459 family: 0.73 in GaPL), justified by the 80.5% of total phenotypic variation explained by PC1, 460 while PC2 accounted only for 12.7%. Despite the lower variability of PC2, this trait could be 461 predicted with a reasonable accuracy of 0.42 in the collection, while in most of the families 462 the accuracy level was above 0.2 (with, and in some cases without TRS optimization). PC2 463 was not predictable in 'GDFj' and 'FjPL', two families with moderate and high transgression 464 on the PC2 axis. The results showed that using PC1 and PC2 as a first tentative to perform a 465 multi-trait prediction was a relevant method to predict fruit texture profiles through an 466 integrative approach. 467 468

Impact of genetic clustering and relatedness on prediction accuracy 469
Having highly related individuals between the TRS and the TS is necessary but not 470 always sufficient for an optimal TRS design; in fact enlarging the TRS with scarcely related 471 individuals can diminish prediction accuracies (Lorenz & Smith, 2015). Moreover, trait 472 variation can be coupled with genetic structure. Several studies have for instance showed the 473 impact of genetic structure on genomic prediction, demonstrating that taking genetic structure 474 into account can improve GS efficiency (Guo et al., 2014;Isidro et al., 2015;Rio et al., 475 2019). Although in apple the genetic structure is known to be weak, with substantial levels of 476 admixture in apple cultivars (Urrestarazu et al., 2016;Vanderzande et al., 2017;Cornille et 477 al., 2019), it could still have a relevant effect on predictions, depending on the population 478 composition and the trait under investigation. Significant genetic structure has been identified, 479 for instance, between dessert and cider apples, which could potentially be correlated with fruit 480 quality traits (Lassois et al., 2016). Through the implementation of the DAPC method, six 481 significant although lowly differentiated genetic clusters were obtained, with families 482 belonging to one or two specific clusters, depending mostly to the assignment of their parental 483 cultivars. While some degree of correlation was apparent between the genetic clustering of 484 individuals and their phenotypic distribution (Fig. S3), the addition of the clustering effect 485 into the prediction model almost systematically degraded the prediction accuracies. Moreover, 486 the TRS optimization based on clustering was the lowest performing among the four methods 487 tested. This could indicate that additive relationship alone already captured the genetic 488 clustering present in our population. One important information given by the clustering 489 patterns was that the 'GaPL' and 'GaPi' families, for which both parents were in the same 490 genetic cluster or in the best represented cluster in the collection (Cluster 5), yielded the best 491 predictions. 492 The genetic parameter having the largest impact on predictions was genetic 493 relatedness, with the two families most related to the collection ('GaPL' and 'GaPi') yielding 494 by far the highest accuracies compared to the remaining Fuji-related families. This 495 observation finds consistency to the fact that genetic relationship is a fundamental parameter 496 in genomic prediction (see e.g. Habier et al., 2010;Clark et al., 2012;Daetwyler et al., 2014). 497 The addition of closely-related individuals from the same family (scenario 2) or from a 498 complete half-sib family (scenario 3) to the collection did not improve the prediction 499 accuracy, except for 'FjPi', for which scenario 2 was the most accurate. This result might 500 indicate that either the collection retains already 'enough' diversity to predict families, or that 501 the excess of unrelated individuals in the collection cannot be corrected by adding related 502 individuals. Thus, scenario 2 and 3 do not seem to effectively improve the TRS. 503 To this end, the gradual increase of the TRS size using a priori information of genetic 504 parameters was used as an alternative optimization strategy. TRS optimization was tested in 505 four different ways, based on a priori information on similarities between individuals. These 506 were represented either by additive relationship or by genetically derived principal 507 components coordinates (Fig. 5, Fig. S7, Table 3, Table S7). The results allowed in all cases 508 to improve predictions tested beforehand with TRS scenarios 1 to 3 with a minimal increase 509 of 0.2 and maximal increase of 0.4, reaching a maximum accuracy of 0.81 ('GaPL', PC1, 510 Table 3). This means that the maximum accuracies were also never reached by employing the 511 entire collection, especially for families with the lowest genetic relatedness to the TRS (i.e. to 512 the collection here). The best prediction accuracy for fruit texture in apple was obtained with 513 the implementation of 50 individuals in the TRS for families less related to the entire TRS and 514 at least 100 accessions for families with a higher genetic relationship (or clustering within the 515 major genetic cluster of the TRS, such as 'GaPL' and 'GaPi' here). These results are 516 consistent with previous findings in barley from Lorenz and Smith (2015), that showed the 517 detrimental effects of adding unrelated individuals to the TS into the TRS, partially 518 contradicting the idea that having at least one related individual in the TRS is sufficient to 519 increase accuracies (Daetwyler et al., 2014). 520 Our results thus provided useful information for the TRS composition, illustrating the 521 complex roles of structure and relatedness in shaping texture variability in apple. 522 523 Towards a simplified assessment of fruit texture for genomic selection 524 The improvement of fruit texture is still limited by the time-consuming and expensive 525 assessment needed for its dissection and the low variation observed in modern elite apple 526 accessions due to the fixation of PG1 (Atkinson et al., 2012;Di Guardo et al., 2017). Thus, 527 even though we demonstrate the feasibility of GS for apple texture, its application will be 528 considered only if predictions are precise enough to perform the costly phenotyping of the 529 TRS. The characterization of texture is a challenging task, as this trait is composed of 530 mechanical and acoustic sub-traits. The analysis of PC1 and PC2 relied on the texture 531 dissection and the measurements of these 12 traits. In particular, FNP, which is the number of 532 mechanical peaks observed in the mechanical profile generated by fruit compression on the 533 texture analyzer, was highly correlated with the group of the acoustic traits related to 534 crispiness. As mechanical traits are easier to measure than acoustic ones, FNP would be in 535 practice the best measurement to choose for assessing crispiness. Since we also obtained high 536 prediction accuracy for FNP (0.63 in collection and maximum of 0.78 in "optimized" family 537 prediction), we propose this sub-trait as the most valuable descriptor for fruit texture, 538 minimizing the effort needed to phenotype such as complex phenotype. Moreover, the 539 predictions presented in this study have been performed with a set of 8,294 SNPs, which is 540 still not dense enough considering the rapid decay of the linkage disequilibrium in apple 541 (Laurens et al. 2018). Although we reached already satisfying accuracies with this amount of 542 SNP, it would be useful to increase the number of markers with the available apple 480K 543 (Bianco et al., 2016) or by using genotyping-by-sequencing methods  to 544 further improve predictions. 545 The use of principal component as synthetic traits resulted to be a valuable multi-trait 546 approach to better predict and understand the texture variability. Here, we investigated in 547 details that fruit crispiness (PC2) in particular is less variable than fruit firmness (PC1). While 548 crispy apples are necessarily firm, the opposite relationship is in fact not validated. Our 549 predictions indicated that fruit firmness in apple can be accurately selected (along PC1), but it 550 needs to be taken into account that an excessive value for this trait can lead to unpleasant 551 quality perception for the consumer. On the other side, crispiness was better predicted with 552 the PC2. Despite the lower variation for crispiness in our population, the selection for this trait 553 resulted to be feasible, although with lower accuracy. To improve the predictions for 554 crispiness we might need to increase the variability for this trait within the TRS. More 555 generally, while the selection on fruit traits has shaped apple domestication, the current 556 cultivated pool relies on a few founders, hence having a narrow genetic basis. Thus, a better 557 targeting of apple texture might necessitate a pre-breeding step incorporating or generating 558 genetic diversity for this trait with the use of mealy cultivars and of wild relatives of Malus 559 domestica (Khan et al., 2014;Peace et al., 2019). 560 561 Supplementary data 562 Table S1. Texture genotypic values and coordinates for PC1 and PC2. 563 Table S2. Additive relationship matrix. 564 Table S3. Assignments of individuals to genetic clusters. 565 Table S4. Pairwise Fst-values between genetic clusters. 566 Table S5. Accuracies obtained in cross-validations within the collection using two models. 567 Table S6. Accuracies obtained in family predictions using two models and three TRS 568 scenarios. 569 Table S7. Accuracies obtained in family predictions with TRS optimization with four 570 methods. 571