The Product neutrality function defining genetic interactions emerges from mechanistic models of cell growth

  1. Lucas Fuentes Valenzuela
  2. Paul Francois
  3. Jan M Skotheim  Is a corresponding author
  1. Department of Biology, Stanford University, United States
  2. Department of Biochemistry and Molecular Medicine, University of Montreal, Canada
  3. Chan Zuckerberg Biohub, United States

eLife Assessment

The paper addresses the question of gene epistasis and asks what is the correct null model for which we should declare no epistasis. By reanalyzing synthetic gene array datasets regarding single and double-knockout yeast mutants, and considering two theoretical models of cell growth, the authors reach the valuable conclusion that the product function is a good null model. While the justification of some assumptions is incomplete, the results have the potential to be of value to the field of gene epistasis.

https://doi.org/10.7554/eLife.105265.3.sa0

Abstract

Genetic analyses, which examine the phenotypic effects of mutations both individually and in combination, have been fundamental to our understanding of cellular functions. Such analyses rely on a neutrality function that predicts the expected phenotype for double mutants based on the phenotypes of the two individual non-interacting mutations. In this study, we examine fitness, the most fundamental cellular phenotype, through an analysis of the extensive colony growth rate data available for budding yeast. Our results confirm that the Product neutrality function describes the colony growth rate, or fitness, of a double mutant as the product of the fitnesses of the individual single mutants. This Product neutrality function performs better than Additive or Minimum neutrality functions, supporting its continued use in genetic interaction studies. Furthermore, we explore the mechanistic origins of this neutrality function by analyzing two theoretical models of cell growth. We perform a computational genetic analysis to show that in both models, the Product neutrality function naturally emerges due to the interdependence of cellular processes that maximize growth rates. Thus, our findings provide mechanistic insight into how the Product neutrality function arises and affirm its utility in predicting genetic interactions affecting cell growth and proliferation.

Introduction

Genetic analysis has been one of the primary methods scientists use to understand how a cell works. One way this is done is through the analysis of genetic interactions, in which the phenotypic effects of mutations are analyzed both as single mutations and then together as double mutations in the same cell. Genes are then said to interact if their combination produces phenotypes that are different from what is predicted from a generic model combining two non-interacting mutations (Phillips, 2008; Beltrao et al., 2010; Costanzo et al., 2019). In other words, a genetic interaction is identified when combining multiple mutations yields something unexpected. Yet, what should we expect when combining mutations in such a complex system as a living cell?

The expected phenotype of a double mutant predicted from the two single mutants’ phenotypes is defined by the neutrality function. In this way, the neutrality function calculates the expected phenotype of a double-mutant strain carrying two non-interacting mutations (Beltrao et al., 2010; Mani et al., 2008). If the double-mutant phenotype deviates significantly from that given by the neutrality function for two specific mutations, they are then said to interact. A lot of care, therefore, needs to be taken in selecting the appropriate neutrality function, which depends on the context and phenotype to be examined (Phillips, 2008). The neutrality function should be defined such that most mutations are categorized as non-interacting. If the neutrality function were not defined this way, the majority of genes would appear to interact, leaving only a few genes with distinct functions. However, decades of cell biological and structural biological analysis have identified specific functions for many genes and their associated proteins. For example, the components of the ribosome or RNA polymerase have the specific task to form these complexes, and metabolic enzymes catalyze specific biochemical reactions. This implies that a judiciously selected neutrality function should predict most double-mutant fitnesses from the single-mutant fitnesses since any two randomly selected genes should be unlikely to interact. Recent technological advances have enabled the screening of the proliferation of single and double genetic mutants at increasingly larger scales in E. coli (Typas et al., 2008; Butland et al., 2008; Babu et al., 2014), fission yeast (Roguev et al., 2008; Dixon et al., 2008), C. elegans (Lehner et al., 2006; Byrne et al., 2007), and human cells (Horlbeck et al., 2018). In budding yeast, Synthetic Genetic Arrays (SGAs) Tong et al., 2001 have generated the largest and most comprehensive such datasets (Baryshnikova et al., 2010b; Costanzo et al., 2010; Costanzo et al., 2016).

The most fundamental phenotype of a cell is its fitness, namely how quickly it grows and proliferates in a given environment. Here, we focus on this property and define fitness as the relative exponential growth rate with respect to that of a wild-type cell. For yeast proliferation, single-mutant fitnesses are usually assumed to combine according to a Product neutrality function, namely that the fitness of the double mutant is the product of the fitnesses of the single mutants (Mani et al., 2008). This neutrality function was shown to better predict double-mutant fitnesses than an Additive neutrality function, where the differences between wild-type and mutant fitnesses were simply added together, and a Minimum neutrality function, where the double-mutant fitness was taken as the lowest fitness of the two single-mutant strains. However, this analysis was based on older, less extensive data, which raises the question of whether this neutrality function remains accurate when the large amount of more recently collected yeast data is also considered. And, if so, then why does the Product neutrality function accurately describe double-mutant fitnesses? In other words, what are the properties of the underlying system controlling cell growth and proliferation that result in a Product neutrality function for mutant fitnesses?

In this paper, we conduct an in-depth analysis of recent yeast double-mutant datasets and show that they support the Product neutrality function. Moreover, we analyze two theoretical models of growth of increasing complexity (Scott et al., 2010; Weiße et al., 2015) and find that the Product neutrality function emerges naturally from both growth models, albeit with small deviations specific to each. Taken together, our work supports the use of the Product neutrality function to model genetic interactions in the regulation of cell growth and gives mechanistic insight into its origin.

Results

High-throughput gene perturbation experiments in budding yeast support a Product neutrality function for double-mutant fitness

To test the general validity of neutrality functions for mutations affecting cell proliferation, we sought to examine the most extensive such dataset. In the SGA dataset, over 20 million single- and double-mutant budding yeast strains were generated. Then, the growth rates of their colonies were measured in SGAs (Tong et al., 2001; Baryshnikova et al., 2010b; Costanzo et al., 2010; Costanzo et al., 2016). Each mutant’s fitness was then defined as this measured growth rate normalized by that of the wild-type strain, enabling a consistent comparison across thousands of genotypes. Next, we use SGA datasets of growth of single- and double-mutant cells for pairs of mutations to test specific neutrality functions.

Here, we followed (Mani et al., 2008) and began our examination using the Product, Additive, and Minimum neutrality functions (see Figure 1A). The Product neutrality function predicts that the fitness of a double mutant is the product of the fitnesses of the two corresponding single mutants, while the Additive neutrality function proposes that the difference between the double-mutant and wild-type fitnesses is the sum of the differences between the two mutant and wild-type fitnesses. The Minimum neutrality function proposes that the fitness of the double mutant is equal to the fitness of the least fit single mutant. Effectively, these neutrality functions express different forms of modularity or independence between cellular processes. The Product model suggests that a mutation’s effect depends on the fitness of the background strain without that mutation, while a mutation’s effect is independent of the fitness of the background strain in an Additive model. The Minimum model suggests that there is some rate-limiting process whose slow time scale dominates the determination of cell growth so that more minor mutations affecting other processes have no additional effect.

High-throughput gene deletion experiments in budding yeast support a Product neutrality function for double-mutant fitness.

(A) Budding yeast mutant fitness is defined as the colony growth rate relative to that of wild-type cells. Schematic illustration of epistasis in growth rate and of different laws proposed in the literature. λ denotes the growth rate, and W the fitness. (B–D) For each double mutant, we plot the residual of the fitness predicted from the indicated model against the fitness of the fittest of the two separate single mutants (maximum single-mutant fitness). Dots indicate the median for 10 equally spaced bins between 0.5 and 1. (E) Box plots for the distributions of the residuals for the three neutrality functions as a function of the maximum single-mutant fitness. A thick line denotes the median, and boxes denote the 25th and 75th percentiles of the distributions. The data plotted here represents a subset of the entire Synthetic Genetic Array (SGA) dataset, corresponding to the Deletion Mutant Array (DMA) at 30°C.Appendix 1—figures 1 and 2 report results for the other subdatasets.

To compare the different neutrality functions with double-mutant fitnesses, we first perform some minor pre-processing of the SGA data (see Methods for details). We note that these data are from combinations of gene deletions, temperature-sensitive alleles, and hypomorphic mutants (see Methods). Moreover, cells were growing quickly on the relatively rich synthetic complete media containing glucose (Baryshnikova et al., 2010a). For these reasons, increasing the growth rate is difficult, and we are in the regime where mutations generally decrease fitness. This contrasts with evolution experiments, where fitness increases very slowly through the gradual accumulation of mutations, which likely exhibit different neutrality functions from those we consider here (Phillips, 2008; Bakerlee et al., 2022; Johnson et al., 2023). Consistent with previous work (Mani et al., 2008), we see that the Product neutrality function better predicts double-mutant fitnesses as a function of the single-mutant fitnesses over a broad range of fitness defects (Figure 1, Appendix 1—figure 1). Indeed, the median residual for the Product neutrality function remains very close to zero even for highly deleterious mutations, while it significantly deviates for the other two. For instance, for a maximum single-mutant fitness of 72%, the median residual is −0.8% for the Product neutrality function, while they are −16.8% and 9.8% for the Minimum and Additive neutrality functions, respectively. For smaller values of the maximum single-mutant fitness, the median residual remains virtually unchanged for the Product neutrality function, while it deviates even further from zero for the other two. However, we note that significant variation around the median residual remains. At a maximum single-mutant fitness of 72%, the interquartile range lies between 8% and 10% depending on the neutrality function considered.

This observation that the Product neutrality function best describes double mutant fitnesses holds for mutations to essential and nonessential genes, and across different temperature conditions (see Methods and Appendix 1—figures 1 and 2). The Minimum neutrality function generally predicts fitnesses that are too high, while the Additive neutrality function generally predicts fitnesses that are too low.

The Product neutrality function describes interactions between genes associated with two distinct biological processes

While the Product neutrality function predicts double-mutant fitnesses better than the other ones we considered, there remains significant variation (residuals) in the data. This suggests that mutations affecting different functional parts of the cell might be following different neutrality functions. To determine whether this is the case, we analyze the distribution of epistasis residuals for pairs of distinct biological processes and their associated genes. We use the Gene Ontology (GO) dataset (Ashburner et al., 2000, Aleksander et al., 2023) and, for each biological process, extract the genes in the SGA dataset that are associated with this process (see Methods for details). In particular, we define inter-process gene pairs as pairs of gene perturbations in two distinct biological processes identified using GO annotations (Figure 2A). Similarly, intra-process pairs are defined as pairs of gene perturbations in the same GO-annotated biological process (Appendix 1—figure 3A). Then, for each pair of GO-defined processes, we compute the residuals for the neutrality functions for all pairs of mutations where one mutation is associated with one process and the second mutation with the other. For each neutrality function and pair of GO processes, we extract the median residual as a function of the largest single-mutant fitness defect (a proxy for mutation severity). This shows that, while imperfect, the Product neutrality function is a good description of typical interactions, while the Additive and Minimum neutrality functions have large, systematic residuals (Figure 2B–D). In general, this result is expected and consistent with the use of this type of genetic analysis to define mutations in genes from different biological processes as not interacting. Moreover, we find that this is not only generally true, but also true for each specific pair of processes that we consider. In other words, we do not find evidence that there are particular pairs of processes whose mutations significantly deviate from the Product neutrality function or more closely follow an Additive or Minimum neutrality function.

The Product neutrality function describes interactions between genes associated with two distinct biological processes.

(A) Schematic illustration of the analysis process. We first select two different Gene Ontology (GO) biological processes and extract the double mutants in the Synthetic Genetic Array (SGA) dataset associated with them. Then, we compute the median residual for each pair of biological processes and each neutrality function. (B–D) Median residual for the Minimum, Product, and Additive neutrality functions as a function of the maximum single-mutant fitness. Each line denotes mutations to a different pair of distinct GO biological processes. The majority of biological process pairs closely follow the Product model.

While mutations associated with different biological processes have fitnesses generally predicted by the Product neutrality function, since they generally do not interact, this may not be the case for mutations associated with the same process. For example, if two mutations break the same protein complex, one would not expect any additional drop in fitness for the double mutant. Consistent with this notion, for gene pairs in the same biological process, we observe more deviations as well as significantly larger residuals (Appendix 1—figure 3B–D). The quantitative comparison of the two types of interactions reveals that large residuals (both positive and negative) are significantly more likely for two mutations categorized as being in the same GO process (Appendix 1—figure 3E, F).

We note that the SGA dataset has already been used to assign a biological process to each gene (Costanzo et al., 2016). For each mutation, a vector of residuals for the Product neutrality function with all other mutations was generated. Then, a Pearson correlation coefficient was calculated for each pair of these vectors. The reasoning was that mutations affecting the same biological process should have similar genetic interaction profiles, which was found to be the case. This then allowed the clustering of groups of correlated mutations, which were named using prior knowledge of many genes in each cluster. That this analysis generally works, that is, the gene clusters have discernible biological meaning, can be viewed as further support of the Product neutrality function.

A bacterial growth model partially supports the Product neutrality function

Having verified empirically that the Product neutrality function is supported by the latest data for cell proliferation, we now turn our attention to its origins. Addressing this question requires some mechanistic model of biosynthesis. However, most mechanistic models of growth apply directly to single cells in rich nutrient conditions, which may not directly apply to the SGA measurements of colony expansion rates. In particular, colony growth has been shown to follow a biphasic pattern (Meunier and Choder, 1999). A first exponential phase is followed by a slower linear phase as the colony expands. Previous modeling and empirical work indicates that this second linear expansion rate reflects the underlying exponential growth of cells in the periphery of the colony (Pirt, 1967; Gray and Kirwan, 1974; Baryshnikova et al., 2010a; Gandhi et al., 2016; Zackrisson et al., 2016; Miller et al., 2022). More precisely, mathematical models show the linear colony-size expansion rate is directly proportional to the square root of the exponential growth rate under non-limiting conditions. Intuitively, this relationship arises because colony growth is dominated by the expansion of the population of cells in an annulus at the colony border that are exposed to rich nutrient conditions. These cells expand at a rate similar to the exponential rate of cells growing in a rich nutrient liquid culture. In contrast, the cells in the interior of the colony experience poor nutrient conditions, grow very slowly, and do not contribute to colony growth.

This intimate relationship between both proliferation rates allows us to explore the origin of the Product neutrality function in mechanistic models of cell growth. Indeed, if colony-based fitnesses follow a Product model, then

WxycWxcWycλxycλWTcλxcλyc(λWTc)2,

where the superscript c indicates colony-based values for the fitness W and the growth rate λ. Taking into account the relationship between single-cell exponential growth rates and colony growth rates, we can write

λcλl,

where the superscript l denotes liquid cultures. Combining these expressions, we obtain

λxylλWTlλxlλylλWTl2WxylWxlWyl.

In other words, from the perspective of the Product neutrality function, fitnesses based on colony expansion rates are equivalent to fitnesses based on single-cell exponential growth rates. The prevalence of the Product neutrality model—both in the SGA data and in previous studies on datasets from liquid cultures (Mani et al., 2008; Jasnos and Korona, 2007; Onge et al., 2007)—encourages the exploration of its origin in mechanistic models of cell growth.

While models of entire cells do exist, these models are complex and computationally intensive (Onge et al., 2007; Karr et al., 2012). This makes probing and extracting explanatory information from these models difficult. We therefore sought to analyze simpler, more tractable, lower-dimensional models of cell growth. Coarse-grained models offer an appealing alternative for probing the fundamental principles of metabolism and growth (Scott et al., 2010; Weiße et al., 2015; Roy et al., 2021; Balakrishnan et al., 2022; Chure and Cremer, 2023; Calabrese et al., 2023). Rather than representing as many reactions as possible, they provide an integrated representation of generic processes in the cell. Their simplicity and low dimensionality make them easy to compare with empirical measurements and to examine for potential explanatory relationships.

The reduced, tractable model of cell growth that we will consider first was developed for E. coli (Scott et al., 2010; Scott and Hwa, 2011; Figure 3A). While this model was developed for E. coli bacteria and validated using data from this organism, there is nothing specific to prokaryotes in the model. Experimental measurements in other organisms suggest that the observations leading to this model, including that the cellular ribosome fraction increases with growth rate, are in fact generic and also seen in the yeast S. cerevisiae (Metzl-Raz et al., 2017; Elsemman et al., 2022; Xia et al., 2022). In its simplest form, the model defines growth as resulting from two sets of processes, metabolic and translational, that interact in a linear pathway. The metabolic sector provides precursors that are then assembled into proteins by the translational sector, and the flux through each sector is determined by the amount of proteins in that sector. For optimal growth, in which no proteins are wasted, the flux through the metabolic sector is equal to the flux through the translational sector so that

λ=κtϕt=κnϕn,
A bacterial growth model partially supports the Product neutrality function.

(A) Schematic of the bacterial growth model by Scott and Hwa. Growth rate is defined by the translation flux, which is itself equal to the metabolic flux. The cell partitions its proteome so as to maximize growth rate. (B) Mutations are modeled such that they affect either of the parameters, separately. Values of κt and κn in the mutant are indicated with primes and are sampled from a uniform distribution from 0 to their value in wild-type cells. indicates the corresponding growth rate. The analytical expression of the double-mutant fitness consists of the Product model with a perturbation. (C–E) For each sampled double mutant, we plot the residual of the fitness predicted from the indicated model against the fitness of the fittest of the two separate single mutants (maximum single-mutant fitness). Dots indicate the median for 10 equally spaced bins between 0.5 and 1. (F) Box plots for the distributions of the residuals for the three models and the model in (C) as a function of the maximum single-mutant fitness. A thick line denotes the median, and boxes denote the upper and lower quartiles of the data. The analytical model in (B), named Scott–Hwa and shown in red, is exact.

where the flux through translation and metabolic sectors is characterized by the parameters κt and κn multiplying the fraction of the proteome devoted to ribosomes and metabolic proteins, ϕt and ϕn, respectively. These fluxes directly determine the growth rate of the cell, λ. The total proteome is fixed and is partitioned into metabolic, translational, and ‘other’ sectors of the cell so that

1=ϕt+ϕn+ϕo

where ϕo is the fraction of the cell devoted to other housekeeping functions. This set of algebraic equations can be solved for the optimal growth rate

λ=(1ϕo)11/κt+1/κn.

We can then model a mutation as a perturbation to these parameters that decreases the growth rate since these are the types of mutations that dominate the budding yeast data (Figure 3B). Single mutants have either κt or κn perturbed, while double mutants have both parameters perturbed. Given that we are concerned with two non-interacting mutations, we do not consider the cases where both mutations affect the same parameters. Indeed, one expects two mutations affecting the same parameter to interact. These combinations are therefore inappropriate to study the emerging neutrality functions from growth models, and we do not consider them in this paper.

We can also consider mutations to ϕo, which could be associated with deleting a gene encoding a protein that is not required for the given growth condition. This would serve to increase the cell growth rate because now a larger fraction of the proteome could be devoted to metabolism and translation. In this case, a mutation to ϕo and another to either κt or κn would combine exactly multiplicatively so that

Wxy=WxWy.

However, in the generally rich media conditions the SGA experiments were done, there is no evidence that any gene deletion causes an increase in cell growth rate so we do not consider this type of mutation further. We therefore ignore the multiplicative factor (1ϕo) and analyze the following expression for growth rate

λ=11/κt+1/κn.

We can then analytically derive a closed-form solution for the double-mutant fitness as a function of the single-mutant fitnesses (Figure 3C; see SI for details). Under this model, which we call Scott–Hwa in reference to the authors of the initial work, we observe that the Product neutrality function fits the mutational analysis of the Scott–Hwa model better than the Additive or Minimum neutrality functions (Figure 3D-G).

We understand the better performance of the Product neutrality function to arise from a type of feedback regulation that ensures that the flux through all sectors is equal. This effectively makes the mutations to the metabolic and translational sector interdependent, despite their a priori independent functions. In this model, these sectors are coupled because the cell is assumed to have a feedback process to optimize growth rate under any perturbation to its parameters. For instance, in the case of a mutation that decreases the metabolic capacity κn, this feedback drives an increase in the fraction of metabolic proteins at the expense of translational proteins such that the growth rate is maximized under those new parameters. If this feedback were absent, then the growth rates of the double mutant would be significantly lower. In that case, the double-mutant fitness actually follows a Minimum neutrality function (see Appendix 1 and Appendix 1—figure 4). Finally, we also note that the Product neutrality function does not accurately predict model fitnesses well for beneficial mutations—that is, mutations that increase growth rate (see Appendix 1—figure 5). This is because the deviations from the Product neutrality function in the mathematical derivation for the double-mutant fitness in Figure 3 can diverge for beneficial mutations (e.g. Wx=Wy=2). When the Scott–Hwa model’s growth-optimizing feedback operates in the context of beneficial mutations, as one process is made more efficient, proteomic resources are allocated to accelerate other processes in the cell. In this way, improving the efficiency of one process will indirectly benefit other processes, leading to compound effects such that double-mutant fitnesses are higher than any of the three models predicts for beneficial mutations.

We note that in this analysis, we do not aim to replicate the statistics of mutations in the SGA dataset, where mutations to either sector could be statistically rarer or more frequent than the other. Instead, we here aim to analyze how mutations to independent parameters governing cell growth combine considering the simplest model with two sectors and their corresponding parameters.

The Product neutrality function accurately predicts fitness for many pairs of parameters in a more complex cell growth model

While the Scott–Hwa model has proven successful for predicting many aspects of bacterial growth, it remains very simple. Therefore, we sought to explore a more complex model that explicitly incorporates more aspects of biosynthesis. Here, we consider the model of Weiße et al., 2015, which incorporates nutrient intake, transcription, competitive binding between mRNAs and ribosomes, and translation, all of which are mediated by associated enzymes and a limiting cellular ‘energy’ (Figure 4A). The Weiße model decomposes cell growth into multiple steps (see supporting material for a full model description). External nutrients are first imported into the cell and then metabolized into a cellular ‘energy’. Both of these steps are catalyzed by associated transport and metabolic enzymes according to Michaelis–Menten kinetics. Transcription and translation are then activated by this generated ‘energy’, also via Michaelis–Menten kinetics. In particular, the model incorporates different transcription rates for ribosomal and non-ribosomal mRNAs. Different mRNAs then compete for free ribosomes to form a ribosome–mRNA complex. This mRNA competition is modeled using mass action kinetics with specified binding and unbinding rates. Four types of proteins are explicitly modeled as a product of translation: transport proteins, metabolic enzymes, ribosomal proteins, and so-called q-proteins which support housekeeping functions much like the ‘other’ proteins in the Scott–Hwa model.

The Product neutrality function accurately predicts fitness for many pairs of parameters in a more complex cell growth model.

(A) Schematic of the growth model from Weiße et al., 2015. This model includes nutrient intake, metabolism, transcription, and translation. (B) Schematic of the mutational analysis. For each pair of parameters α and β, mutations are modeled such that they affect either of the parameters, separately. Then, the median residual is computed for each neutrality function and they are subsequently reported for every pair of parameters considered. For each parameter pair, we report the mean deviation of the simulated double mutants from (C) the Product neutrality function and (D) the analytical expression of the double-mutant fitness under the Scott–Hwa model. Only parameter pairs corresponding to two different biological processes are considered. Those corresponding to the same process are grayed out. Parameter pairs involving translation (gmax,Kp) are the ones described best by the Scott–Hwa model, while the others are better described by the Product neutrality function.

To perform a mutational analysis of the Weiße model, we first identified the parameters in the model that can reasonably be expected to change through a gene perturbation (see supporting material for details). For instance, we assume that some parameters, such as the maximum nutrient import rate, could be impacted by a gene perturbation, while other parameters, such as the average gene length in the genome, could not. This led to the identification of 9 easily interpreted parameters whose mutation could negatively impact the cell growth rate and correspond to 28 parameter pairs associated with different biological processes. We then performed a similar analysis as we did for the Scott–Hwa model. Namely, we constructed mutants for each pair of parameters by rescaling the values of the original parameters by a number randomly sampled between 0 and 1 (see Methods for details). We then analyzed the statistics of the epistasis coefficients for the neutrality functions that we considered so far (Figure 4B).

Our mutational analysis of the Weiße model revealed several striking observations. First, the Product neutrality function is generally better than the Additive or Minimal neutrality functions at describing the mutational results. However, we observe a range of responses and can identify two key subpopulations of parameter pairs. The subset of parameter pairs involving protein translation follows the Scott–Hwa model very closely (Figure 4D). This is not entirely surprising, as the Weiße model is an extension of the Scott–Hwa model and incorporates a similar global feedback optimizing cell growth and a competition for resources. On the other hand, another subset of parameter pairs follows the Product neutrality function even more closely. These parameter pairs involve mutations to the other sectors, including the transport, metabolism, and transcription sectors (Figure 4C). This raises the question of why some parameter pairs more closely follow the Product neutrality function than others.

Nonlinear kinetics drive deviations from the Product neutrality function in the Weiße model for cell growth

To address the question as to what drives deviations from the Product neutrality function in our genetic analysis of the Weiße model, we took an analytical approach. We examined the dependence of the growth rate on the parameter pairs exhibiting small deviations from the Product neutrality function. To do this, we first extracted a closed form expression that models the growth rate λ and its dependence on two mutated parameters α and β (Figure 5A; see supporting material). While only approximate, this derivation represents the data appropriately for members of this subset of parameters (see Appendix 1—figure 7). There are two striking features in this derivation. First, the growth rate λ has an explicit dependence on the two parameters α and β, when α and β are selected from the subset of parameters governing metabolism and transport sectors. Second, the amplitude of deviation from the Product neutrality function is governed by an additional parameter, γ, which is the inverse of the Michaelis–Menten constant giving the transcription rate as a function of the cellular ‘energy’. Thus, when γ is small, transcription is a less efficient process that is then linearly related to the cellular ‘energy’ available. When γ is larger, transcription is saturated and performed at a rate unrelated to the available ‘energy’. As we decrease γ, we observe that the Product neutrality function is a better and better approximation (Figure 5C). Importantly, the intuition provided by the analytical approximation extends to multiple pairs of parameters (see Appendix 1 and Appendix 1—figure 8). Taken together, our analysis of the Weiße model shows how the Product neutrality function naturally arises for many different parameter pairs and how deviations from it can be driven by nonlinear effects, such as those that can emerge from Michaelis–Menten kinetics.

Nonlinear kinetics drive deviations from the Product neutrality function in the Weiße model.

(A) A subset of parameter pairs we analyzed follows the Product model very closely. We derived an analytical approximation of the growth rate and the double-mutant fitness for these pairs and found that the deviation from the product law is governed by nonlinear kinetics. (B) In the case of the parameter pair (vt,ns), we show that the deviation from the Product model is driven by the Michaelis–Menten constant θx associated with transcription (see Supporting Information). (C) Tuning the value of γ impacts how good of an approximation the Product neutrality function is for this and other parameter pairs (see text). This analysis validates the analytical approximation and highlights how nonlinear kinetics, in this case Michaelis–Menten kinetics, can drive deviations from the Product neutrality function.

Discussion

Cell growth and proliferation are fundamental to cell biology and have been subject to extensive genetic analysis aiming to understand the underlying regulatory network. Such genetic analysis often aims to identify interactions between mutations through the combination of individual mutations in a double-mutant cell. If the double-mutant proliferates at an unexpected rate, the mutated genes are considered to interact. This, of course, raises the question as to what is the expected rate of proliferation for a cell containing both mutations given the proliferation rate of a cell containing only one of the individual mutations. By analyzing a high-throughput dataset of interactions of gene perturbations in budding yeast (Costanzo et al., 2010; Costanzo et al., 2016), we found that single-mutant fitnesses tend to combine multiplicatively, consistent with earlier work (Mani et al., 2008).

After establishing that the fitness of a double mutant is expected to be approximately the product of the fitness of the individual mutants, namely, the Product neutrality function, we sought to determine if this was also a feature of models of cell growth. If so, then what underlying mechanisms present in these models give rise to this Product neutrality function? Our analysis complements previous, more abstract theoretical attempts at understanding the origin of the Product neutrality function that are not based on any specific model of cell growth (Chiu et al., 2012). Indeed, we found that the Product neutrality function best fits the simulated double-mutant fitnesses despite deviations that depend on the specific parameter pairs and the particular model considered.

That the Product neutrality function fits the budding yeast data and cell growth models better than the Additive and Minimum neutrality functions has important implications for the underlying network controlling cell growth and proliferation. On the one hand, the Minimum neutrality function implies that the double-mutant fitness is set by the most deleterious mutation so that the process this gene is involved in becomes rate limiting for cell growth. Clearly, cell growth does not operate this way, likely because the underlying processes are interconnected. Mutations impairing protein translation impact synthesis of all the proteins in the cell so that other processes, like transcription or surface transport, are also affected. As the cell readjusts its machinery to ensure optimal use of resources, interconnected processes are impacted through a redistribution of cellular resources. On the other hand, the Additive neutrality function implies that a mutation has the same absolute effect on the proliferation rate regardless of the presence of another mutation. This is also clearly not the case as the Additive neutrality function consistently predicts fitnesses below those observed in the data. In the models, this is due in part to growth-supporting feedback that reapportions the proteome. In reality, this may reflect the presence of the general stress response which supports cells in response to genetic or environmental perturbations limiting their growth rate (Gasch et al., 2000). In this way, the Product neutrality function is a reasonable intermediate model between Minimum and Additive that incorporates—albeit approximately—effects such as growth-optimizing feedback and is consistent with the phenomenon of diminishing returns or increasing cost epistasis (Reddy and Desai, 2021). Moreover, our theoretical analysis gives insight into the mechanistic underpinning of the Product neutrality function. In our analyses, a product of the single-mutant growth rates naturally emerges in the analytical treatment of both theoretical models that we consider, albeit with deviation terms that depend on the specific model and simplifying assumptions.

Taken together, our work here constitutes a first step toward understanding the structure of interactions inherent in cell growth models. While we focused on coarse-grained models for their simplicity and mechanistic interpretability, they might be too simple to effectively model large double-mutant datasets and the resulting double-mutant fitness distributions. For instance, it is not possible to differentiate between multiple types of growth rate perturbations impacting the same sector, as they would all be modeled through a limited number of parameters (Metzl-Raz et al., 2017). We therefore expect the combination of high-throughput genetic data with the analysis of larger-scale models, for instance based on Flux Balance Analysis, Metabolic Control Analysis, or whole-cell modeling, to lead to important complementary insights regarding the regulation of cell growth and proliferation (Karr et al., 2012; Oftadeh et al., 2021; Segrè et al., 2005; He et al., 2010; Orth et al., 2010; Kacser and Burns, 1973; Szathmáry, 1993; Dykhuizen et al., 1987; de Vienne et al., 2023; Kryazhimskiy, 2021) We also believe that theoretical exploration of fitness landscapes will shed light on the underlying structure of growth and metabolism networks (Reddy and Desai, 2021; Guo et al., 2019; Boffi et al., 2023). In addition to larger-scale models, we see the refinement of the measurement of cell growth rates as a path forward to a better understanding of its regulation (MacLean, 2010). While we showed here that the Product neutrality function fits the data well for deleterious mutations, we anticipate that there are significant and meaningful deviations that are currently obscured by experimental noise. Similarly, large-scale measurements of the impacts of beneficial mutations will be instrumental in testing the validity of the Product neutrality function in this other regime. From our modeling efforts, we anticipate that such measurements could give important insights into the underlying genetic network regulating growth and proliferation.

Methods

Analysis of the SGA dataset

The complete SGA dataset was accessed on the cell map webpage. In the SGA genetic interaction dataset, a set of query mutant strains is crossed to an ordered array of mutants.

There are two sets of query mutants. The first one consists of a mix of nonessential deletion mutant strains and of temperature-sensitive alleles of essential genes. The second one is a set of mutants carrying hypomorphic, Decreased Abundance by mRNA Perturbation (DAmP) alleles of essential genes.

There are also two types of arrays. The Deletion Mutant Array (DMA) denotes deletions to a set of nonessential genes, while the Temperature Sensitive Array (TSA) contains a mix of essential and nonessential genes.

Both sets of query mutants are crossed to either type of array, at two different temperature conditions, namely 26 and 30°C. The analysis of Figure 1 reports the analysis of the first set of query mutants crossed to the DMA at 30°C. In Appendix 1—figure 1, we report the same analysis for the first set of query mutants for the other array–temperature combinations. In Appendix 1—figure 2, we report the analysis for the DAmP set of query mutants in the different array–temperature combinations.

For each subdataset—that is, each combination of query mutants, array, and temperature condition—the data is processed in the following steps. First, only deleterious mutations are kept. That is, we remove mutants having a fitness larger than 1. Second, we eliminate mutants where the Additive model predicts a negative fitness, that is, such that

Wx+Wy<1,

because in this case the prediction under the Additive neutrality function is negative (Wx+Wy1<0). While we could have analyzed these datapoints with the other neutrality functions, we sought to analyze all neutrality functions on the same consistent dataset.

For the sake of clarity, the scatter plots in Figure 1 do not reproduce the entire dataset. Instead, the dataset is binned in 10 bins along the x-axis, and 500 values are sampled at random in that bin. However, the median lines on top of the scatter plots (e.g. Figure 1B–D) as well as the box plots (e.g. Figure 1, Appendix 1—figures 1 and 2) apply to the entire dataset.

Analysis of the GO biological processes

The analysis of GO biological processes is based on the Uniprot database. In this dataset, genes are associated with a series of GO biological processes. To analyze the behavior of the fitness of double mutants associated with different biological processes, we first selected the set of biological processes that were represented by a large enough number of single mutants in the SGA dataset. Arbitrarily, this limit was set at 50. This led to a limited number of 47 biological processes (and a maximum total of 1081 pairs), which we report in the section ‘Analysis of GO biological processes’. Naturally, as genes are potentially associated with multiple GO biological processes, this analysis sometimes leads to pairs of processes with genes in common. In this case, we discarded the pair so that we only consider biological process pairs that do not have any genes in common. This results in 685 pairs of biological processes that do not share any genes.

Mutational analysis of growth models

A mutation is modeled as a perturbation of a parameter that decreases the growth rate. For a given parameter α, we model a perturbation as α=θα, where θ is a random variable uniformly distributed in [0,1]. When estimating the impact of the parameter γ in the Weiße model (see section ‘Nonlinear kinetics drive deviations from the Product neutrality function in the Weiße model for cell growth’), we perform a mutational analysis as described above for different values of the parameter γ and collect the median residual.

Appendix 1

A.1 Analysis of GO biological processes

The 47 biological processes that result from the analysis (with their GO code) are given in Appendix 1—table 1.

A.2 Scott–Hwa model

A.2.1 Derivation of the double-mutant fitness

Here, we derive the expression for the double-mutant fitness in the Scott–Hwa model. First of all, the wild-type growth rate λWT is written as

(1) λWT=κtκnκt+κn.

As we consider only deleterious mutations, we can model a perturbation to parameters as

(2) κn=(1ε)κn,
(3) κt=(1ε)κt,

where δ,ε[0,1]. For simplicity, we will denote by the x subscript mutations affecting κt and the y subscript mutations affecting κn. We note the order does not matter. Therefore, the different fitnesses W=λ/λWT are given by

(4) Wx=1ε1εκt/(κt+κn),
(5) Wy=1δ1δκn/(κt+κn),
(6) Wxy=(1ε)(1δ)(κt+κn)(1ε)κt+(1δ)κn.

From there, we can calculate that

(7) WxWy=(1ε)(1δ)(κt+κn)(1ε)κt+(1δ)κn+εδλWT.

This expression is similar to Equation 6. Indeed, we can rewrite it as

(8) Wxy=WxWy(1+εδλWT(1ε)κt+(1δ)κn).

We see that the double-mutant fitness Wxy consists of the product of single-mutant fitnesses and a deviation. We will now express this deviation as a function of single-mutant fitnesses only, showing that the expression does not depend on the value of the parameters κt and κn.

Rearranging Equations 4 and 5, we have

(9) 1Wx=εκn(1ε)κt+κn,
(10) 1Wy=δκt(1δ)κn+κt,
(11) (1Wx)(1Wy)=1κt+κnεδλWT(1εκt/(κt+κn))(1δκn/(κt+κn)).

which can be rearranged to see that

(12) εδλWT(1ε)κt+(1δ)κn=(1Wx)(1Wy)(1εκt/(κt+κn))(1δκn/(κt+κn))(1ε)κt+(1δ)κn
(13) =(1Wx)(1Wy)1εκt/(κt+κn)δκn/(κt+κn)+εδλWT/(κt+κn)(1ε)κt+(1δ)κn
(14) =(1Wx)(1Wy)(1ε)κt+(1δ)κn+εδλWT(1ε)κt+(1δ)κn
(15) =(1Wx)(1Wy)(1+εδλWT(1ε)κt+(1δ)κn).

We see a recursive relationship appear, such that

(16) εδλWT(1ε)κt+(1δ)κn=(1Wx)(1Wy)(1+(1Wx)(1Wy)(1+))
(17) =(1Wx)(1Wy)+((1Wx)(1Wy))2+((1Wx)(1Wy))3+

Therefore, we have from Equation 8 that

(18) Wxy=WxWyk=0((1Wx)(1Wy))k
(19) =WxWy11(1Wx)(1Wy)
(20) =WxWy(1+(1Wx)(1Wy)1(1Wx)(1Wy)),

from the property of geometric series and noting that (1Wx)(1Wy)1.

A.2.2 Scott–Hwa model with no feedback

To understand the impact of growth rate optimization in the Scott–Hwa model, we consider an alternative model formulation where this feedback is absent. We assume that we still have two different processes (metabolism and translation) and that the flux through both of them has to be equal. However, we do not incorporate the assumption that the proteome has to be partitioned between these two processes. Instead, each process is assumed to have a given amount of proteins or enzymes available at its disposal. This formulation implies that

(21) λ=κtϕt=κnϕn,
(22) 0ϕtϕtmax,
(23) 0ϕnϕnmax,

where κt and κn are the capacities of the translation and metabolic sector, respectively, as in the Scott–Hwa model, and where ϕt and ϕn are the normalized concentrations of proteins or enzymes.

If we assume that the cell maximizes its growth rate, we have that

(24) λmaxnf=min(κtϕtmax,κnϕnmax),

where the nf superscript stands for no feedback.

Clearly, in this simplified model, one sector will be limiting either because it is not efficient enough, that is, κi is too small, or because it does not have enough resources to deploy to accelerate the process, that is, ϕimax is too small. In this case, we see that the growth rate λ corresponds to the minimum of the potential maximal fluxes in either sector.

Therefore, modeling mutations as affecting either κt or κn but not the protein fractions, we find that

Wx=min(κtϕtmax,κnϕnmax)λ0,Wy=min(κtϕtmax,κnϕnmax)λ0,Wxy=min(κtϕtmax,κnϕnmax)λ0

For the sake of argument, assume that κtϕtmaxκnϕnmax. Therefore, λ0=κtϕtmax and Wx=κtϕtmaxλ0=κtκt, and Wy=min(1,κnϕnmaxκtϕtmax).

We can then write that

Wxy=min(Wx,κnϕnmaxκtϕtmax).

If κnϕnmaxκtϕtmax>1, then Wxy=Wx1. However, if κnϕnmaxκtϕtmax1, then Wy=κnϕnmaxκtϕtmax. Therefore, we have that

Wxy=min(Wx,Wy),

that is, in the model with no feedback laid out above, double-mutant fitnesses are actually the minimum of the single-mutant fitnesses. Simulation results are reported in Appendix 1—figure 4.

We note that this absence of feedback results in significantly smaller growth rates. In the Scott–Hwa model, we have

λSH=κtκnκt+κn.

This implies that ϕt=κnκt+κn and ϕn=κtκt+κn.

Let us now consider, in the context of the model in this section, ϕtmax=κnκt+κn and ϕnmax=κtκt+κn. Upon mutation of κt or κn, we have

(25) λnfλSH=min(κtκnκt+κn,κnκtκt+κn)κtκnκt+κn
(26) =min(κn(κt+κn)κn(κt+κn),κt(κt+κn)κt(κt+κn))
(27) =κt+κnκt+κnmin(κnκn,κtκt)

Assume, without loss of generality, that κtκtκnκn, that is κnκnκtκt. Then, we have

κt+κnκt+κnκt+κtκnκtκt+κn=κtκt,

which implies that λnf/λSH1. Therefore, in the absence of feedback, the growth rate is smaller than in the Scott–Hwa model, which is consistent with the interpretation of growth-rate optimization.

A.3 Weiße model

A.3.1 Original model

The original model consists of a system of equations describing the synthesis, degradation, and reaction of a series of molecules and proteins. In particular, external nutrients s are converted into internal nutrients si. Those nutrients are then converted into a generic cellular energy a that enables transcription and translation. Indeed, mRNAs m are transcribed and then bind to ribosomes to form an mRNA–ribosome complex c that governs the rate of protein translation. There are four main types of mRNAs and complexes in the model, indicated by a subscript. Each is associated with different processes: transport, metabolism, ribosomal, and q-proteins. This last type of protein denotes housekeeping proteins whose concentration remains approximately constant.

(28) si˙=νimp(et,s)νcat(em,si)λsi
(29) a˙=nsνcat(em,si)x{r,t,m,q}nxνx(cx,a)λa
(30) r˙=νr(cr,a)λr+x{r,t,m,q}(νx(cx,a)kbrmx+kucx)
(31) et˙=νt(ct,a)λet
(32) em˙=νm(cm,a)λem
(33) q˙=νq(cq,a)λq
(34) m˙x=ωx(a)(λ+dm)mx+νx(cx,a)kbrmx+kucx
(35) c˙x=λcx+kbrmxkucxνx(cx,a)

All of these reactions are modulated by parameters reported in Appendix 1—table 2.

The rates governing the system of equations are assumed to follow Michaelis–Menten kinetics as follows

(36) νimp(et,s)=etvtsKt+s
(37) νcat(em,si)=emvmsiKm+si
(38) νx(cx,a)cxγ(a)nx
(39) γ(a):=γmaxaKγ+a
(40) ωx(a)=ωxaθx+a

Finally, the growth rate corresponds to total mass of proteins being synthesized at steady state, that is

(41) λ=γ(a)xcxM.

A.3.2 Isolation of parameters

This model has a total of 21 parameters whose values were set (either from the literature or estimated) in the original model (Table S2 in Weiße et al., 2015; reproduced here under Appendix 1—table 2).

Among them, some represent quantities or parameters that would not reasonably change under a mutation. Therefore, we do not consider the following parameters in our mutation analysis: external nutrients s, mRNA degradation rate dm, ribosome length nr, length of non-ribosomal proteins nx, q-autoinhibition Hill coefficient hq, the mRNA–ribosome binding/unbinding rates kb and ku, the total cell mass M.

To determine whether a numerical mutational analysis with the remaining 13 parameters could be considered, we computed the impact of a change of these parameters on the growth rate λ (Appendix 1—figure 6). Upon mutation, most parameters do impact the growth rate negatively. However, two parameters vm,Km do not seem to have an impact on the growth rate. We attribute this to the metabolic sector not being limiting in this model, with those parameter values. We, therefore, exclude those two parameters from further analysis. We also exclude the parameters associated with the q-proteins Kq and ωq as these are associated mostly with housekeeping functions.

We are therefore left with nine parameters modeling the impact of four different processes. In total, this results in 28 potential combinations of two parameters associated with different processes.

A.3.3 Simplification

To facilitate analytical treatment, we consider a simplified model with only two types of proteins, instead of four: a general protein p and ribosomes r. We also eliminate the equations corresponding to the competitive binding of mRNAs for ribosomes. Practically, this assumption means that every mRNA binds immediately to a ribosome. While this assumption might not be true in all regimes, we find that this simplifies the analytical treatment and still enables us to understand something concrete about the model.

The simplified model we consider is

(42) s˙i=νimp(p,s)νcat(p,si)λsi
(43) a˙=nsνcat(p,si)xnxνx(cp,a)λa
(44) r˙=νr(cr,a)λr+xνx(cx,a)
(45) p˙=νp(p,a)λp
(46) c˙r=ωr(a)νr(cr,a)λcr
(47) c˙p=ωp(a)νp(cp,a)λcp

If we assume the dilution fluxes to be negligible in comparison to other fluxes in Equations 42–47, we get that

(48) νimp(p,s)νcat(p,si)pvtsKt+spvmsiKm+si
(49) nsνcat(p,si)xnxνx(cp,a)=(cp+cr)γ(a)
(50) ωp(a)νp(cp,a)
(51) ωr(a)νr(cr,a)

The dilution flux needs to be nonzero for p, as its governing equation only contains two terms. We have therefore that

(52) νp(p,a)=λppωpaθx+a1λ

As the growth rate λ=γ(a)xcxM, we have

(53) λ=nsνcatM1MnsvtpsKt+s

Injecting Equation 52 into Equation 53, we have

(54) λ1MnsvtwpsKt+saθx+a

From the above expression, we see that terms belonging to different sectors combine as a product. However, the energy a is not a parameter, but rather a variable in the model. Its value also directly depends on the values of other parameters.

From Equation 49, we have that

(55) λ1M(npωpaθx+a+nrωraθr+a)

If θx<a<θr, we can approximate that with

(56) λ1M(npωp+nrωraθr).

If aθx,θr, we have

(57) λ1M(npωpθx+nrωrθr)a.

In both cases, we can assume that there is a linear relationship between the growth rate λ and the energy vector a, such that

(58) abλ+c.

We can inject that in the a/(θx+a) term in Equation 54, which gives

(59) aθx+aλ+c/b(θx+c)/b+λ

Assuming that c/bλ, we have

(60) aθx+aλK+λ,

with K=(θx+c)/b.

Injecting the above in Equation 54, we get that,

(61) λ1MnsvtwpsKt+sλK+λ.

Therefore, K acts as a scale of growth rate where deviations from the Product model can appear. If Kλ, then

(62) λ1MnsvtwpsKt+s.

On the other hand, if Kλ, then

(63) λ1MnsvtwpsKt+s.

In both of these cases, the growth rate λ behaves as the product of multiple parameters. However, if λK, and if we denote two parameters in the expression by α, β without loss of generality, we have

(64) λ2(K+λ)αβλ.

Either λ=0 (which we will not consider), or

λ(K+λ)=αβλ2+Kλαβ=0λ=K+K2+4αβ2=K2(1+1+4αβK2)

We can write, for small x, that

1+x=1+x2x28+

Assuming that 4αβK2 is small and keeping terms up to second order, we therefore have that

λ=K2x2(1x4), with x=4αβK2λαβ(1γαβ), with γ=1/K2

A.4 Double-mutant fitness

From there, we can compute an approximation to the double-mutant fitness. Let us assume that λ=αβ(1γαβ) as derived above. To simplify notation, we will rescale the parameters to their WT value, that is

α=α/αWT,β=β/βWT,γ¯=γαWTβWT.

Under this rescaling, we see that

λWT=αWTβWT(1γ¯),Wx=α(1γ¯α)1γ¯,Wy=β(1γ¯β)1γ¯,Wxy=αβ1γ¯αβ(1γ¯).

Noting that

WxWy=αβ(1γ¯α)(1γ¯β)(1γ¯)2

we can show that

Wxy=WxWyαβ(1γ¯)2γ¯(1α)(1β),=WxWy(1γ¯1α1γ¯α1β1γ¯β).

If γ is small, then Wx=α(1γ¯α)1γ¯α and Wy=β(1γ¯β)1γ¯β, and therefore we can write that

WxyWxWy(1γ¯1Wx1γ¯Wx1Wy1γ¯Wy).
Appendix 1—figure 1
Double-mutant fitnesses are best described by the Product neutrality function in the Synthetic Genetic Array (SGA) dataset.

Box plots for the distributions of the residuals for the three neutrality functions as a function of the maximum single-mutant fitness. Each plot corresponds to a different subset of the SGA dataset. Namely, they correspond to the first set of query mutants (see Methods) crossed to different types of mutant arrays in different temperature conditions (A) Deletion Mutant Array at 26°C. (B) Temperature Sensitive Array at 26°C. (C) Temperature Sensitive Array at 30°C. Thick line denotes the median, and boxes denote the 25th and 75th percentiles of the distributions. The Product neutrality function models the data consistently better than the others.

Appendix 1—figure 2
Double-mutant fitnesses are best described by the Product neutrality function in the Synthetic Genetic Array (SGA) dataset.

Box plots for the distributions of the residuals for the three neutrality functions as a function of the maximum single-mutant fitness. Each plot corresponds to a different subset of the SGA dataset. Namely, they correspond to the second set of query mutants (DAmP, see Methods) crossed to different types of mutant arrays in different temperature conditions (A) Deletion Mutant Array at 30°C. (B) Temperature Sensitive Array at 26°C. (C) Temperature Sensitive Array at 30°C. Thick line denotes the median, and boxes denote the 25th and 75th percentiles of the distributions. The Product neutrality function models the data consistently better than the others.

Appendix 1—figure 3
Larger deviations from the Product neutrality function characterize gene pairs affecting the same GO biological process.

(A) Schematic illustration of the analysis process for double mutants where both mutations affect the same GO biological process. We first select two different GO biological processes and extract the double mutants in the Synthetic Genetic Array (SGA) dataset associated with them. Then, we compute the median residual for each pair of biological processes and each neutrality function. (B–D) Median residual for the Minimum, Product, and Additive neutrality functions as a function of the maximum single-mutant fitness. Each line denotes mutations to a single GO biological process. We see larger deviations from the Product model than in Figure 2. (E) Histogram of the SGA dataset after extracting pairs affecting either two different (inter) or the same (intra) GO biological process. Large residuals are much more likely when both mutations affect the same GO biological process.

Appendix 1—figure 4
The Scott–Hwa model with no feedback follows a Minimum neutrality function.

Box plots for the distributions of the residuals for the model of Scott–Hwa model with no feedback as a function of the maximum single-mutant fitness. A thick line denotes the median, and boxes denote the upper and lower quartiles of the data. The absence of feedback due to resource competition in the model of Scott–Hwa model with no feedback results in a Minimum neutrality function.

Appendix 1—figure 5
Large deviations from the Product neutrality function characterize beneficial mutations.

Box plots for the distributions of the residuals for the different neutrality functions for beneficial mutations as a function of the minimum single-mutant fitness. Thick lines denote the median, and boxes denote the upper and lower quartiles of the data. Deviations from the Product neutrality function derived in Derivation of the double-mutant fitness can be unbounded.

Appendix 1—figure 6
Eleven parameters exhibit negative impact on growth rate upon mutation in the Weiße model.

Among the initial 21 parameters, 13 were kept as candidates for a mutational analysis. Two of them vm,Km, associated with the metabolic sector, do not have any impact on growth rate upon mutation, likely because this sector is not limiting for growth in that parameter range. The others have a negative impact upon mutation.

Appendix 1—figure 7
Deviations from the Product neutrality function in the Weiße model are captured by the γ approximation.

Box plots for the distributions of the residuals for all models considered in this paper, for two example parameter pairs (A: ωq,ns, B: Kt,ωe). A thick line denotes the median, and boxes denote the upper and lower quartiles of the data. The Gamma model, in purple, denotes the derivation in Figure 5A. It captures the small deviations from the Product neutrality function.

Appendix 1—figure 8
Tuning γ impacts how good an approximation the Product neutrality function is for multiple parameter pairs in the Weiße model.

The analysis of Figure 5C, illustrating the impact of γ on a single parameter pair (ns,vt) is here extended to other parameter pairs to demonstrate the validity of the mechanistic interpretation. (A) vt,ωr; (B) ωe,ns; (C).θr,vt In all cases, we see that decreasing γ results in better alignment with the Product model.

Appendix 1—table 1
GO biological processes used in the analysis.
GO biological process nameGO biological identifier
DNA integrationGO:0015074
DNA recombinationGO:0006310
DNA repairGO:0006281
Ascospore formationGO:0030437
Cell cycleGO:0007049
Cell divisionGO:0051301
Cell wall organizationGO:0071555
Cellular response to DNA damage stimulusGO:0006974
Cellular response to oxidative stressGO:0034599
Chromatin remodelingGO:0006338
Chromatin silencing at telomereGO:0006348
Chromosome segregationGO:0007059
Cytoplasmic translationGO:0002181
EndocytosisGO:0006897
Endoplasmic reticulum to Golgi vesicle-mediated transportGO:0006888
Fungal-type cell wall organizationGO:0031505
Intracellular protein transportGO:0006886
Intracellular signal transductionGO:0035556
mRNA splicing, via spliceosomeGO:0000398
MacroautophagyGO:0016236
Maturation of SSU-rRNA from tricistronic rRNA transcriptGO:0000462
Meiotic cell cycleGO:0051321
Mitochondrial translationGO:0032543
Negative regulation of transcription by RNA polymerase IIGO:0000122
Positive regulation of transcription by RNA polymerase IIGO:0045944
Proteasome-mediated ubiquitin-dependent protein catabolic processGO:0043161
Protein foldingGO:0006457
Protein import into nucleusGO:0006606
Protein phosphorylationGO:0006468
Protein targeting to vacuoleGO:0006623
Protein transportGO:0015031
Protein ubiquitinationGO:0016567
Pseudohyphal growthGO:0007124
rRNA methylationGO:0031167
rRNA processingGO:0006364
Reciprocal meiotic recombinationGO:0007131
Regulation of transcription by RNA polymerase IIGO:0006357
Regulation of transcription, DNA-templatedGO:0006355
Ribosomal large subunit biogenesisGO:0042273
Sporulation resulting in formation of a cellular sporeGO:0030435
Transcription by RNA polymerase IIGO:0006366
Transcription elongation from RNA polymerase II promoterGO:0006368
Translational terminationGO:0006415
Transmembrane transportGO:0055085
Transposition, RNA-mediatedGO:0032197
Ubiquitin-dependent protein catabolic processGO:0006511
Vesicle-mediated transportGO:0016192
Appendix 1—table 2
Model parameters from Weiße et al., 2015, obtained either from the literature or from parameter optimization.
DescriptionDefault valueUnit
sExternal nutrient104[molecs]
dmmRNA-degradation rate0.1[min−1]
nsNutrient efficiency0.5None
nrRibosome length7459[aa/molecs]
nx,x{t,m,q}Length of non-ribosomal proteins300[aa/molecs]
γmaxMax. transl. elongation rate1260[aa/min molecs]
KγTransl. elongation threshold7[molecs/cell]
vtMax. nutrient import rate726[min−1]
KtNutrient import threshold1000[molecs]
vmMax. enzymatic rate5800[min−1]
KmEnzymatic threshold1000[molecs/cell]
wrMax. ribosome transcription rate930[molecs/min cell]
we=wt=wmMax. enzyme transcription rate4.14[molecs/min cell]
wqMax. q-transcription rate948.93[molecs/min cell]
θrRibosome transcription threshold426.87[molecs/cell]
θnrNon-ribosomal transcription threshold4.38[molecs/cell]
Kqq-Autoinhibition threshold152,219[molecs/cell]
hqq-Autoinhibition Hill coeff.4None
kbmRNA–ribosome binding rate1[cell/min molecs]
kumRNA–ribosome unbinding rate1[min−1]
MTotal cell mass108[aa]

Data availability

The current manuscript is a computational study, so no data have been generated for this manuscript. Modeling code is uploaded to the Skotheimlab Github repository https://github.com/skotheimlab/GrowthModels (copy archived at Fuentes Valenzuela, 2025).

The following previously published data sets were used

References

    1. Aleksander SA
    2. Balhoff J
    3. Carbon S
    4. Cherry JM
    5. Drabkin HJ
    6. Ebert D
    7. Feuermann M
    8. Gaudet P
    9. Harris NL
    10. Hill DP
    11. Lee R
    12. Mi H
    13. Moxon S
    14. Mungall CJ
    15. Muruganugan A
    16. Mushayahama T
    17. Sternberg PW
    18. Thomas PD
    19. Van Auken K
    20. Ramsey J
    21. Siegele DA
    22. Chisholm RL
    23. Fey P
    24. Aspromonte MC
    25. Nugnes MV
    26. Quaglia F
    27. Tosatto S
    28. Giglio M
    29. Nadendla S
    30. Antonazzo G
    31. Attrill H
    32. Dos Santos G
    33. Marygold S
    34. Strelets V
    35. Tabone CJ
    36. Thurmond J
    37. Zhou P
    38. Ahmed SH
    39. Asanitthong P
    40. Luna Buitrago D
    41. Erdol MN
    42. Gage MC
    43. Ali Kadhum M
    44. Li KYC
    45. Long M
    46. Michalak A
    47. Pesala A
    48. Pritazahra A
    49. Saverimuttu SCC
    50. Su R
    51. Thurlow KE
    52. Lovering RC
    53. Logie C
    54. Oliferenko S
    55. Blake J
    56. Christie K
    57. Corbani L
    58. Dolan ME
    59. Drabkin HJ
    60. Hill DP
    61. Ni L
    62. Sitnikov D
    63. Smith C
    64. Cuzick A
    65. Seager J
    66. Cooper L
    67. Elser J
    68. Jaiswal P
    69. Gupta P
    70. Jaiswal P
    71. Naithani S
    72. Lera-Ramirez M
    73. Rutherford K
    74. Wood V
    75. De Pons JL
    76. Dwinell MR
    77. Hayman GT
    78. Kaldunski ML
    79. Kwitek AE
    80. Laulederkind SJF
    81. Tutaj MA
    82. Vedi M
    83. Wang S-J
    84. D’Eustachio P
    85. Aimo L
    86. Axelsen K
    87. Bridge A
    88. Hyka-Nouspikel N
    89. Morgat A
    90. Aleksander SA
    91. Cherry JM
    92. Engel SR
    93. Karra K
    94. Miyasato SR
    95. Nash RS
    96. Skrzypek MS
    97. Weng S
    98. Wong ED
    99. Bakker E
    100. Berardini TZ
    101. Reiser L
    102. Auchincloss A
    103. Axelsen K
    104. Argoud-Puy G
    105. Blatter M-C
    106. Boutet E
    107. Breuza L
    108. Bridge A
    109. Casals-Casas C
    110. Coudert E
    111. Estreicher A
    112. Livia Famiglietti M
    113. Feuermann M
    114. Gos A
    115. Gruaz-Gumowski N
    116. Hulo C
    117. Hyka-Nouspikel N
    118. Jungo F
    119. Le Mercier P
    120. Lieberherr D
    121. Masson P
    122. Morgat A
    123. Pedruzzi I
    124. Pourcel L
    125. Poux S
    126. Rivoire C
    127. Sundaram S
    128. Bateman A
    129. Bowler-Barnett E
    130. Bye-A-Jee H
    131. Denny P
    132. Ignatchenko A
    133. Ishtiaq R
    134. Lock A
    135. Lussi Y
    136. Magrane M
    137. Martin MJ
    138. Orchard S
    139. Raposo P
    140. Speretta E
    141. Tyagi N
    142. Warner K
    143. Zaru R
    144. Diehl AD
    145. Lee R
    146. Chan J
    147. Diamantakis S
    148. Raciti D
    149. Zarowiecki M
    150. Fisher M
    151. James-Zorn C
    152. Ponferrada V
    153. Zorn A
    154. Ramachandran S
    155. Ruzicka L
    156. Westerfield M
    157. Gene Ontology Consortium
    (2023) The gene ontology knowledgebase in 2023
    GENETICS 224:iyad031.
    https://doi.org/10.1093/genetics/iyad031
    1. Kacser H
    2. Burns JA
    (1973)
    The control of flux
    Symposia of the Society for Experimental Biology 27:65–104.

Article and author information

Author details

  1. Lucas Fuentes Valenzuela

    Department of Biology, Stanford University, Stanford, United States
    Contribution
    Conceptualization, Formal analysis, Investigation, Writing – original draft, Writing – review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0009-0002-7403-3537
  2. Paul Francois

    Department of Biochemistry and Molecular Medicine, University of Montreal, Montreal, Canada
    Contribution
    Conceptualization, Formal analysis, Supervision, Investigation, Writing – review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-2223-839X
  3. Jan M Skotheim

    1. Department of Biology, Stanford University, Stanford, United States
    2. Chan Zuckerberg Biohub, San Francisco, United States
    Contribution
    Conceptualization, Resources, Supervision, Funding acquisition, Investigation, Writing – original draft, Project administration, Writing – review and editing
    For correspondence
    skotheim@stanford.edu
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0001-8420-6820

Funding

National Institutes of Health (GM134858)

  • Jan M Skotheim

The funders had no role in study design, data collection, and interpretation, or the decision to submit the work for publication.

Version history

  1. Sent for peer review:
  2. Preprint posted:
  3. Reviewed Preprint version 1:
  4. Reviewed Preprint version 2:
  5. Version of Record published:

Cite all versions

You can cite all versions using the DOI https://doi.org/10.7554/eLife.105265. This DOI represents all versions, and will always resolve to the latest one.

Copyright

© 2025, Fuentes Valenzuela et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 1,249
    views
  • 14
    downloads
  • 0
    citations

Views, downloads and citations are aggregated across all versions of this paper published by eLife.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Lucas Fuentes Valenzuela
  2. Paul Francois
  3. Jan M Skotheim
(2025)
The Product neutrality function defining genetic interactions emerges from mechanistic models of cell growth
eLife 14:RP105265.
https://doi.org/10.7554/eLife.105265.3

Share this article

https://doi.org/10.7554/eLife.105265