^{1}

^{1}

^{2}

^{1}

Conceived and designed the experiments: HB. Performed the experiments: RR. Analyzed the data: RR HB. Contributed reagents/materials/analysis tools: RR MEA HB. Wrote the paper: RR HB. Provided additional advice for statistical analysis: MEA.

Ever since the case of the missing heritability was highlighted some years ago, scientists have been investigating various possible explanations for the issue. However, none of these explanations include non-chromosomal genetic information. Here we describe explicitly how chromosomal and non-chromosomal modifiers collectively influence the heritability of a trait, in this case, the growth rate of yeast. Our results show that the non-chromosomal contribution can be large, adding another dimension to the estimation of heritability. We also discovered, combining the strength of LASSO with model selection, that the interaction of chromosomal and non-chromosomal information is essential in describing phenotypes.

Genome-wide association studies (GWAS) have contributed to the identification of many human loci associated with a wide range of complex traits, such as height, intelligence or diseases such as obesity, type 2 diabetes and age-related macular degeneration. However, GWAS do not explain the whole story of the observed heritability of these traits [

Interestingly, non-chromosomal genetic information has not yet been taken into account when estimating heritability, although there is evidence for effects on the phenotype arising from cytoplasmic elements in many different organisms. For instance, Cadwell et al. [

Recently, our collaborators from MIT have designed an experiment using yeast where both sources of information, chromosomal and non-chromosomal were controlled in order to observe the phenotype of a single chromosomal polymorphism in the presence of different cytoplasmic elements [

Unfortunately, their statistical analysis had some limitations. Firstly, they split the data into training and test sub-samples (1/10 of the data held-out for testing) in order to conduct ten-fold cross validation. This appears to be the wrong approach, given that the sample sizes (see supplementary ^{2}), here used as metric for recovered heritability, for three different models that consider (i) only the chromosomal mutation, (ii) the effects of the chromosomal mutation and non-chromosomal information, and (iii) both the effect of the chromosomal mutation and non-chromosomal information as well as their interaction. They inferred—exclusively from the gain in ^{2}—that non-chromosomal information and interaction effects substantially contribute to the heritability. However, it is well known that ^{2} values may increase with increasing number of explanatory variables, and, hence, can not be exclusively applied to meaningfully compare models with different number of variables.

In this study, we applied more sophisticated statistical means and models, such as the adjusted coefficient of determination (

Previous work studied two inherited non-chromosomal modifiers. First, the presence and absence of the endogenous dsRNA yeast “killer” virus that is transmitted by mitosis and meiosis [^{S288c}) or Sigma mitochondria ([rho+]^{Sigma}). The mitochondria differ strongly in their genomes, with about 2–3 SNPs per kilobase and ten times more insertions and deletions compared to the chromosomal genome. The studied model chromosomal modifiers were gene deletions in the yeast strains S288c and Sigma [

We analysed the effects of non-chromosomal modifiers on growth phenotypes of chromosomal variants using raw data from Edwards et al. [_{ij} as the colony size of controlled genotype i in replicate j. We applied a variance-stabilising Box-Cox transform on the _{ij} values with exponent of 0.25 to obtain normalised values _{ij}, as this parameter choice gave the least-correlated means and variances in the analysed sets (according to Edwards et al. [_{0} + _{1}_{1} + _{0} + _{1}_{1} + _{2}_{2} + _{0} + _{1}_{1} + _{2}_{2} + _{3}_{1}_{2} + _{1}, …, _{n})^{t} as the response vector, _{1}, …, _{n})^{t} ∼ ^{2}_{n}) as the noise vector, _{j} as the _{1}, …, _{p})^{t} as the vector of parameters of interest to be estimated. Each _{j}, _{j} and the response

The ^{2}, because it is capable of handling the inflation of ^{2}, when comparing different models.

In

Three linear models with different complexity are applied to measure the fraction of phenotypic variance. The first model (

Additional to the

908 | 81 | 11 | |

0 | 0 | 1000 | |

0 | 0 | 1000 | |

60 | 121 | 819 | |

3 | 13 | 984 | |

32 | 34 | 934 | |

35 | 31 | 934 | |

606 | 212 | 182 | |

18 | 4 | 978 | |

1 | 18 | 981 |

Each row shows the frequency of selected (according to the MSE) linear models for a gene deletion experiment within 1000 modelling repeats. The interaction model is, aside from the control

Despite these results indicating the importance of non-chromosomal information, we applied LASSO along BIC to verify our findings.

We analysed the effects of the three predictors (_{1}, _{2} and _{1}_{2}) using LASSO alongside BIC. We take advantage of the fact that model selection criteria extracts from a set of candidate models (in our case the

As in the previous approach, we performed 1000 LASSO regressions with BIC model selection for each gene deletion experiment. We then studied the complexity of the BIC-selected models that best describe the given data. In

429 | 544 | 27 | 0 | |

0 | 0 | 451 | 549 | |

0 | 11 | 442 | 547 | |

0 | 618 | 155 | 227 | |

0 | 0 | 992 | 8 | |

0 | 13 | 576 | 411 | |

0 | 2 | 983 | 15 | |

883 | 106 | 10 | 1 | |

0 | 0 | 891 | 109 | |

0 | 0 | 961 | 39 |

Model selection using BIC was applied in each of the 1000 modelling repeats for all gene deletion experiments. Except for the control

Furthermore, we analysed the frequency of predictors representing chromosomal (_{1}) and non-chromosomal effects (_{2}) as well as their interaction (_{1}_{2}) within the 1000 procedures (see _{1} and _{1}_{2} are selected in most cases. Regarding the gene deletion _{1}_{2}) is chosen more often as the other predictors.

_{1} | _{2} | _{1}_{2} | |
---|---|---|---|

518 | 13 | 67 | |

991 | 558 | 1000 | |

907 | 629 | 1000 | |

379 | 237 | 993 | |

1000 | 8 | 1000 | |

997 | 414 | 987 | |

1000 | 5 | 998 | |

89 | 26 | 14 | |

1000 | 109 | 1000 | |

1000 | 39 | 1000 |

Aside from the _{1}) as well as the interaction effect (_{1}_{2}) are selected in almost all modelling processes. An exception is the _{1}_{2}) is most frequently selected.

As for the previous method, these results highlighted the importance of the interaction between chromosomal and non-chromosomal modifiers. However, we discovered that LASSO selected the main or interaction effects arbitrarily (not only for gene deletion experiment

In order to use LASSO for hierarchical interactions, we used the R software package

Interestingly, we identified that in most cases, aside from the experiments _{1}, _{2} and _{1}_{2}) are required to best describe the data (see

Our analyses revealed that the phenotype of a chromosomal mutation may be affected by non-chromosomal elements such as mitochondria and viral state. We also showed that the introduction of non-chromosomal information and its interaction with chromosomal elements considerably enhanced the fraction of explained phenotypic variance of a trait, which is ensured by conserving the chromosomal contribution and the environment, whilst changing the non-chromosomal effects. Previous studies [

However, it is well known that the coefficient of determination, here used as metric for recovered heritability, is prone to an increase when adding more variables to the statistical model. Hence, we could not exclude that the gain in

Due to this fact, we applied model selection criteria BIC to investigate the importance of both chromosomal and non-chromosomal information, as well as their interaction in describing colony size data. BIC highlighted that not only chromosomal mutations, but also the interactions between chromosomal and non-chromosomal elements, such as mitochondria and dsRNA virus are important. The mitochondrial background plays a crucial role in the ^{Sigma}) compared to cells with a S288c mitochondria ([kil-k] [rho+]^{S288c}). With BIC model selection, the _{1}_{2}) is sufficient to describe the data for this gene deletion experiment in two out of three cases.

Furthermore, we examined LASSO with hierarchical interactions. This approach, unlike ordinary LASSO, prevents the inclusion of interaction effects unless the main effects are also included. We identified for most cases, that not only non-chromosomal, but also its interaction with chromosomal effects, is essential to best describe the colony size data. For the

In summary, all applied statistical methods point to the fact that non-chromosomal modifiers, and the interaction effects of chromosomal and non-chromosomal elements, account for a substantial fraction of phenotypic variance of growth rates in yeast.

The raw data measurements and the analysis script (R code) can be found in the supplementary section.

We consider the standard multiple linear regression model with _{1}, …, _{n})^{t} is the response vector, _{1}, …, _{n})^{t} ∼ ^{2}_{n}) is the noise vector; for _{j} represents the _{1}, …, _{p})^{t} is the vector of parameters of interest to be estimated; each _{j}, _{j} and the response

Define ^{2} or the percentage of variance explained is defined as

Since RSS always decreases as more variables are added to the model, ^{2} always increases as more variables are added. For a least squares model with

Maximising

If we have enough observations, we can divide our data set into two parts: a training set of size _{train}, on which the model is fitted, and a test set of size _{test} for evaluation of the performance. A measure of prediction performance commonly used is the _{i, test} and

The LASSO coefficients,

The LASSO technique penalises the regression coefficients using an _{1} norm. It shrinks the coefficients towards zero. In addition, the _{1} penalty has the effect of forcing some of the coefficient to be exactly equal to zero when the tuning parameter

For the least squares model with

Bien el al. [

The additive part is called the ^{p} and Θ ∈ ℝ^{p × p}, where Θ = Θ^{t} and _{jj} = 0 for _{1} = ∑_{j ≠ k}∣_{jk}∣, _{i} is the observed value of _{i}, _{0} is the intercept, and

(PDF)

Click here for additional data file.

(PDF)

Click here for additional data file.

(PDF)

Click here for additional data file.

(R)

Click here for additional data file.

(CSV)

Click here for additional data file.

The yeast experimental design descriptions and the spore growth data were kindly provided by our colleagues David Gifford and Gerald Fink from the Whitehead and Broad Institutes in Cambridge, USA.

Furthermore, we thank Mohammed Dehbi from QBRI and Christopher Leonard from QScience, Qatar Foundation, for improving the quality of the manuscript.