(high level) and Col-0 (low level). The coefficient for this metabolic marker in the C24-heterosis models is negative, which indicates higher heterosis potential of RILs with lower salicylic acid levels. This suggests that the salicylic acid level must not be too high in at least one parent of a well performing testcross. It could be interpreted in such a way that too high alertness of the plant with respect to pathogens (as indicated by too high salicylic acid levels) draws too much of the metabolism into defense associated pathways, thus leading to suboptimal growth. Qualitative prediction of heterosis is possible using the combined genetic-metabolic marker set detected by our quantitative approach As described in the Results section, combined genetic-metabolic models yield a correlation of approximately 0.5 between observed and in LOOV predicted response (cf. Table 1). Although this is a highly significant correlation, as proven by the fact that none of 1000 random permutations of the observed response showed such a high correlation to the predicted response, the respective prediction models cannot be used for reliable quantitative predictions. However, an enrichment of crosses displaying a higher likelihood for strong heterosis would probably be helpful. To evaluate whether this approach might be feasible, we trained a linear discriminant model for C24-heterosis. This model discrim- inates the 359 RILs into lines with high and low heterotic potential using our minimal genetic-metabolic predictor set. In LOOV 233 RILs (i.e., 65%) were classified correctly by that model. Incorporating all 181 metabolic markers and all 110 genetic markers resulted in a discriminant model that classified only 49% of the RILs correctly, which corresponds to a random classifier. Thus the two-step approach described here, i.e. identifying in a quantitative model the most important set of genetic and metabolic markers and subsequently applying those in a qualitative model, might be an interesting path also with respect to practical applications [32]. In addition to the incorporation of gene interactions by metabolic markers, as shown in this study, we will also consider nonlinear relationships of the most important metabolic markers identified in this work. This seems to be a promising step, since on the one hand, the highly predictive metabolic markers show distributions deviating from normal ones. On the other hand, for linear models normal distributions of the predictors are required. Consequently, the consideration of some nonlinear terms in the discrimination model, which account for the special distributions of the highly predictive markers, will probably improve the classification success. Materials and Methods Data and software Data was derived from a RIL population [23] within the framework of Arabidopsis thaliana heterosis experiments that included the accessions C24 and Col-0 [17]. We performed a comparative study on a set of 359 RILs and the two corresponding sets of testcrosses. All variables of interest were available for each of these 1077 lines, and they are described below. The mean frequency of C24 alleles and Col-0 alleles per RIL was 48% and 50%, respectively (cf. Figure S3). In average, 2% of the markers were heterozygous. We investigated two sets of predictor variables, both, indepen- dently and in combination. The first set was the metabolic profile of the RILs consisting of 181 metabolic markers. Their levels were measured via gas-chromatography/mass-spectrometry at 15 days Heterosis Prediction Markers after sowing, and the data was normalized as described elsewhere [19]. Close to one half of the metabolic compounds were identified and 98 substances remained unknown but could at least partly be classified. The second predictor set was composed of 110 DNA- markers distributed along the whole RIL genome [18]. The available genetic marker information had been coded in the following way. Depending on whether the particular DNA-marker featured C24 or Col-0 it was quantified with 1 or 0 while markers with a heterozygous specificity were translated into 0.5. Single missing values have been replaced by the average of the direct neighboured markers’ values in that RIL. If the missing value was at the margin, it just has been replaced by its single neighbour. Markers with more than ten per cent missing values had been deleted in advance. Our response variables of interest were C24-heterosis and Col- heterosis. These were determined by the shoot biomass of the particular RIL, the shoot biomass of the corresponding testcross C24 x RIL (or Col-0 x RIL) and that of the background line C24 (or Col-0). Shoot biomass was measured as dry weight at 15 days after sowing. The determination of the shoot biomass values (cf. Table S6) via mixed model has been described elsewhere [17]. The absolute mid-parent heterosis is defined as the difference between the hybrid’s shoot biomass and the mean shoot biomass of its parents. Dividing the absolute mid-parent heterosis by the mean shoot biomass of the parents, results in the relative mid- parent heterosis, which we focus on in this analysis (cf. Figure 1). In our case, the first parent is C24 (or Col-0) and the second parent is one of the 359 RILs. Our computations were conducted on the free software environment R for statistical computing and graphics [33] in the version 2.5.1. Especially, we used the R package pls [34] in the version 2.0-1. Latent variables For training the particular prediction models, we used the PLS approach [20,35], which is a multiple linear regression method that works on latent variables. It projects the original predictor data into a latent lower dimensional space and thus maximizes the covariance with the response variable. The corresponding R function within the pls package is called plsr [34]. We realise the optimization according to the number of latent variables via LOOV (i.e. N-fold cross-validation with N = 359) of the Pearson correlation between the observed and predicted response. This is achieved in R by setting the parameter validation of the function plsr to ‘‘LOO’’. The correlation, achieved with the optimal number of latent variables, is referred to, in this work, as the predictive power of the considered set of predictors for the response of interest. The R function to accomplish these computations is called plsr. The rational of using PLS instead of ordinary least squares regression is the potential existence of latent variables, on which the particular response depends, but which are not available, because they cannot be measured. Fortunately, one latent variable often regulates not only the response of interest but also several other variables, which are in fact measured in experiments. Since they are all influenced by the one latent variable, many of them correlate strongly. Therefore we can down-weight some of the measured variables. The remaining variables could be considered as biomarkers for the particular response. That means there is not necessarily a causal dependence, but the biomarkers act as an indicator for the response. Since we consider very complex response effects, such as heterosis, probably several latent variables are involved. Thus, it is possible that two different latent variables influence the same measurable variable. For example, one of them up- and the other one down-regulates the measurable variable. PLoS ONE | www.plosone.org 9 April 2009 | Volume 4 | Issue 4 | e5220