Supplementary Notes for:

L. Kuepfer et al., Ensemble modeling as a novel concept to analyze cell signaling dynamics.

Detailed description of methods and additional analysis results.1.Model structures

Each model in the original model ensemble consisted of the core model (CM) and one of the elementary extensions specifying additional control mechanisms. The composition of models in the original and the decoupled ensemble, respectively, is summarized in Table S1. The corresponding biochemical reactions are compiled in Table S2. To decouple the core model, we employed the reaction sets listed in Table S3, whereas the third ensemble was generated by adding elementary extension 11 (Table S1).

Table S1 Characteristics of elementary extensions representing putative biochemical features of the TOR pathway that were used to expand the original core model. Except for the core model, the model extensions were arranged in the order of decreasing overall model quality (increasing normalized training error) for the original model ensemble. Usage of extensions in the original ensemble (OE) and the decoupled ensemble (DE) is indicated by the rightmost columns. The model structure of the decoupled core model (CMD) a priori included elementary extensions 7 and 14. Elementary extensions 5, 10, 12 and 13 were not further considered because the underlying mechanisms were only weakly supported by experimental evidence.

Code / Description / OE / DE
CM / Core model / X / X
1 / Tip41p has two phosphorylation sites / X / X
2 / Tap42p~p-Pph21/22p forms an anti-phosphatase protecting phosphoproteins / X / X
3 / Complex formation of Tap42p~p and Tip41p / X / X
4 / Complex formation of Tap42p~p and Tip41p~p / X / X
5 / Pph21/22p is phosphorylated by Tor1/2p and dephosphorylated by PP2A1/2 / X
6 / Tap42p~p-Pph21/22p acts as a phosphatase / X / X
7 / Specific catalytic constants for dephosphorylation of Tip41p~p by PP2A1/2 / X
8 / Tap42p~p-Sit4p forms an anti-phosphatase that protects phosphorylated proteins / X / X
9 / Tap42p~p-Sit4p is a phosphatase / X / X
10 / Sit4p is phosphorylated by Tor1/2p and dephosphorylated by PP2As / X
11 / Tap42p has two phosphorylation sites / X
12 / PP2A1/2 form with Sit4p and Pph21/22p bound to Tap42~p and dephosphorylate it / X
13 / Tap42p~p bound to Sit4p or Pph21/22p can be dephosphorylated by PP2A1/2 / X
14 / Specific constants for dephosphorylation of Tap42p~p by PP2A1 / PP2A2 / X
15 / Monomeric Sit4p is an active phosphatase for Tip41p~p / X / X
16 / Complex formation of Tap42p and Tip41p~p / X / X
17 / Complex formation of Tap42p and Pph21/22p / X / X
18 / Complex formation of Tap42p and Sit4p / X / X
19 / Tip41p can form a complex with Sap proteins / X

Table S2Biochemical reaction networks for the core model and 18 elementary extensions in the original model ensemble (1-18), and elementary extension 19 (including Sap proteins) that was only employed for the ensemble of decoupled models.

Core model

[Tap42p]+[Tor1/2p] [Tap42p•Tor1/2p]

[Tap42p•Tor1/2p] [Tap42p~p]+[Tor1/2p]

[Tap42p~p]+[Pph21/22p] [Tap42p~p •Pph21/22p]

[Tap42p~p]+[PP2A1] [Tap42p~p •PP2A1]

[Tap42p~p•PP2A1] [Tap42p]+[PP2A1]

[Cdc55p/Tpd3p]+[Pph21/22p] [PP2A1]

[Tap42p~p]+[Sit4p] [Tap42p~p•Sit4p]

[Tap42p~p]+[PP2A2] [Tap42p~p•PP2A2]

[Tap42p~p•PP2A2] [Tap42p]+[PP2A2]

[Sit4p]+[Sap] [PP2A2]

[Tip41p]+[Tor1/2p] [Tip41p•Tor1/2p]

[Tip41p•Tor1/2p] [Tip41p~p]+[Tor1/2p]

[Tip41p~p]+[PP2A1] [Tip41p~p•PP2A1]

[Tip41p~p •PP2A1] [Tip41p]+[PP2A1]

[Tip41p~p]+[PP2A2] [Tip41p~p •PP2A2]

[Tip41p~p •PP2A2] [Tip41p]+[PP2A2]

[Tip41p]+[Tap42p] [Tip41p•Tap42p]

[Fpr1]+[Rapamycin] [Fpr1•Rapamycin]

[Tor1/2p]+[Fpr1p•Rapamycin] [Tor1/2p•Fpr1p•Rapamycin]

1.Tip41p has two phosphorylation sites

[Tip41p~p]+[Tor1/2p] [Tip41p~p•Tor1/2p]

[Tip41p~p•Tor1/2p] [Tip41p~p~p]+[Tor1/2p]

[Tip41p~p~p]+[PP2A1] [Tip41p~p~p •PP2A1]

[Tip41p~p~p •PP2A1] [Tip41p~p]+[PP2A1]

[Tip41p~p~p]+[PP2A2] [Tip41p~p~p •PP2A2]

[Tip41p~p~p •PP2A2] [Tip41p~p]+[PP2A2]

2.Tap42p~p•Pph21/22p forms an anti-phosphatase that protectsphosphorylated proteins from dephosphorylation

[Tap42p~p]+[Pph21/22p] [Tap42p~p •Pph21/22p]

[Tip41p~p]+[Tap42p~p •Pph21/22p] [Tip41p~p •Tap42p~p •Pph21/22p]

[Tap42p~p]+[Tap42p~p •Pph21/22p] [Tap42p~p •Tap42p~p •Pph21/22p]

3. Complex formation of Tap42p~p and Tip41p

[Tip41p]+[Tap42p~p] [Tip41p•Tap42p~p]

4. Complex formation of Tap42p~p and Tip41p~p

[Tip41p~p]+[Tap42p~p] [Tip41p~p •Tap42p~p]

5. Pph21/22p is a phosphoprotein that is phosphorylated by Tor1/2pand dephosphorylated by PP2A1 or PP2A2

[Tor1/2p]+[Pph21/22p] [Tor1/2p•Pph21/22p]

[Tor1/2p•Pph21/22p] [Tor1/2p]+[Pph21/22p~p]

[Pph21/22p~p]+[PP2A1] [Pph21/22p~p •PP2A1]

[Pph21/22p~p •PP2A1] [PP2A1]+[Pph21/22p]

[Pph21/22p~p]+[PP2A2] [Pph21/22p~p •PP2A2]

[Pph21/22p~p•PP2A2] [PP2A2]+[Pph21/22p]

[Cdc55p/Tpd3p]+[Pph21/22p~p] [PP2A1]

6. Tap42p~p•Pph21/22p is a phosphatase that binds tophosphorylated target proteins and dephosphorylates them

[Tap42p~p]+[Pph21/22p] [Tap42p~p •Pph21/22p]

[Tip41p~p]+[Tap42p~p •Pph21/22p][Tip41p~p•Tap42p~p •Pph21/22p]

[Tip41p~p•Tap42p~p•Pph21/22p] [Tip41p]+[Tap42p~p •Pph21/22p]

[Tap42p~p]+[Tap42p~p•Pph21/22p] [Tap42p~p •Tap42p~p •Pph21/22p]

[Tap42p~p •Tap42p~p •Pph21/22p [Tap42p]+[Tap42p~p •Pph21/22p]

7. Specific catalytic constants for the dephosphorylationof Tip41p~p by PP2A1 and PP2A2, respectively

[Tip41p~p •PP2A1] [Tip41p]+[PP2A1]

[Tip41p~p •PP2A2] [Tip41p]+[PP2A2]

8. Tap42p~p•Sit4p is an anti-phosphatase that protectsphosphorylated proteins from dephosphorylation

[Tap42p~p]+[Sit4p] [Tap42p~p •Sit4p]

[Tap42p~p]+[Tap42p~p •Sit4p] [Tap42p~p •Tap42p~p •Sit4p]

[Tip41p~p]+[Tap42p~p •Sit4p] [Tip41p~p •Tap42p~p •Sit4p]

9. Tap42p~p•Sit4p forms a phosphatase that binds tophosphorylated target proteins and dephosphorylates them

[Tap42p~p]+[Sit4p] [Tap42p~p •Sit4p]

[Tap42p~p]+[Tap42p~p •Sit4p] [Tap42p~p •Tap42p~p •Sit4p]

[Tap42p~p •Tap42p~p •Sit4p] [Tap42p]+[Tap42p~p •Sit4p]

[Tip41p~p]+[Tap42p~p •Sit4p] [Tip41p~p •Tap42p~p •Sit4p]

[Tip41p~p•Tap42p~p •Sit4p] [Tip41p]+[Tap42p~p •Sit4p]

10. Sit4p is a phosphoprotein that is phosphorylated by Tor1/2pand dephosphorylated by PP2A1 or PP2A2, respectively

[Tor1/2p]+[Sit4p] [Tor1/2p•Sit4p]

[Tor1/2p•Sit4p] [Tor1/2p]+[Sit4p~p]

[Sit4p~p]+[PP2A1] [Sit4p~p •PP2A1]

[Sit4p~p •PP2A1] [PP2A1]+[Sit4p]

[Sit4p~p]+[PP2A2] [Sit4p~p •PP2A2]

[Sit4p~p •PP2A2] [PP2A2]+[Sit4p]

[Sit4p~p]+[Sap] [PP2A2]

11. Tap42p has two phosphorylation sites

[Tap42p~p]+[Tor1/2p] [Tap42p~p •Tor1/2p]

[Tap42p~p •Tor1/2p] [Tap42p~p~p]+[Tor1/2p]

[Tap42p~p~p]+[PP2A1] [Tap42p~p~p •PP2A1]

[Tap42p~p~p •PP2A1] [Tap42p~p]+[PP2A1]

[Tap42p~p~p]+[PP2A2] [Tap42p~p~p •PP2A2]

[Tap42p~p~p •PP2A2] [Tap42p~p]+[PP2A2]

[Tap42p~p~p]+[Sit4p] [Tap42p~p~p•Sit4p]

[Tap42p~p ~p]+[Pph21/22p] [Tap42p~p~p •Pph21/22p]

12. Cdc55p/Tpd3p can bind to Pph21/22p and Sapp to Sit4p,while being bound to Tap42p~p and dephosphorylate it

[Tap42p~p •Pph21/22p]+[Cdc55p/Tpd3p] [Tap42p~p •PP2A1]

[Tap42p~p •PP2A1] [Tap42p]+[PP2A1]

[Tap42p~p •Sit4p]+[Sap] [Tap42p~p •PP2A2]

[Tap42p~p •PP2A2] [Tap42p]+[PP2A2]

13. Tap42p~p can be dephosphorylated by PP2A1 or PP2A2 whilebeing bound to Sit4p and Pph21/22p

[Tap42p~p •Pph21/22p]+[PP2A1] [Tap42p~p •Pph21/22p•PP2A1]

[Tap42p~p •Pph21/22p•PP2A1] [Tap42p•Pph21/22p]+[PP2A1]

[Tap42p~p •Pph21/22p]+[PP2A2] [Tap42p~p •Pph21/22p•PP2A2]

[Tap42p~p •Pph21/22p•PP2A2] [Tap42p•Pph21/22p]+[PP2A2]

[Tap42~p •Sit4p]+[PP2A1] [Tap42p~p •Sit4p•PP2A1]

[Tap42p~p •Sit4p•PP2A1] [Tap42p•Sit4p]+[PP2A1]

[Tap42p~p •Sit4p]+[PP2A2] [Tap42p~p •Sit4p•PP2A2]

[Tap42p~p •Sit4p•PP2A2] [Tap42p•Sit4p]+[PP2A2]

[Tap42p]+[Pph21/22p] [Tap42p•Pph21/22p]

[Tap42p]+[Sit4p] [Tap42p•Sit4p]

14. Specific catalytic constants for the dephosphorylationof Tap42p~p by PP2A1 and PP2A2, respectively

[Tap42p~p •PP2A1] [Tap42p]+[PP2A1]

[Tap42p~p •PP2A2] [Tap42p]+[PP2A2]

15. Monomeric Sit4p is an active phosphatase for Tip41p~p

[Tip41p~p]+[Sit4p] Tip41p~p •Sit4p]

[Tip41p~p •Sit4p] [Tip41p]+[Sit4p]

16. Complex formation of Tap42p and Tip41p~p

[Tip41p~p]+[Tap42p] [Tip41p~p •Tap42p]

17. Complex formation of Tap42p and Pph21/22p

[Tap42p]+[Pph21/22p] [Tap42p•Pph21/22p]

18. Complex formation of Tap42p and Sit4p

[Tap42p]+[Sit4p] [Tap42p•Sit4p]

19. Tip41p can form a complex with Sap

[Tip41p]+[Sap] [Tip41p •Sap]

[Tip41p~p]+[Sap] [Tip41p~p•Sap]

Table S3 Decoupled reactions from the core model. The corresponding reactions in the original core model (Table S2) were replaced by the reactions listed below.

[Tap42p~p•PP2A1] [Tap42p]+[PP2A1]

[Tap42p~p•PP2A2] [Tap42p]+[PP2A2]

[Tip41p•Tor1/2p] [Tip41p~p]+[Tor1/2p]

[Tip41p~p •PP2A1] [Tip41p]+[PP2A1]

[Tip41p~p •PP2A2] [Tip41p]+[PP2A2]

[Tap42p~p]+[Sit4p] [Tap42p~p•Sit4p]

[Sit4p]+[Sap] [PP2A2]

[Tap42p~p]+[PP2A2] [Tap42p~p•PP2A2]

2.Experimental data used for model development

Previously published data with S. cerevisiae strains W303-1A 1, 2 and JK9-3da/TB50a 3 were considered. All experiments were done in YPD 1, 3 and YPD medium that has been depleted of inorganic phosphate2during the exponential growth phase. For each published data set, a weighting factor M was introduced that reflects the expected accuracy of the experimental method used. Western blot data for Tap42p~p dephosphorylation2 (Table S4) or Tap42p-Sit4p and Tap42p-Tip41p complex formation2, 3 (Table S5, S6) were generally assumed to be ±20% accurate. Relative complex concentrations1, 4 (Table S7) were assumed to be accurate within ±50%. For results that were too low to be properly discriminated during parametrization (e.g. <2% of the overall amount of Pph21/22p associated with Tap42p), the upper bounds were additionally relaxed by allowing a deviation of up to ±100% (Table S7).

The published Tip41p~p dephosphorylation3 and Tap42p-Pph21/22p complex dissociation data1 were not quantitative. As an approximation, the same time course as for Tap42p-Sit4p-protein complex dissociation upon the addition of rapamycin was assumed for the Tap42p-Pph21/22p complex. Relative to the reference state, the Tap42p-Pph21/22p concentration is reduced to 40% at 30 and 60 min after the addition of 109 nM rapamycin to wild type S. cerevisiae 1. Since these percentages were gross estimations, we employed an accuracy of ±100%. The same accuracy was chosen for the dephosphorylation of Tip41p~p, where the level of phosphorylated Tip41p~p was assumed to be 10% of the overall concentration upon addition of 109 nM rapamycin 3.

Absolute protein concentrations 5, 6 were considered to be accurate within ±25% (Table S8). Since the Tor2p concentration was not determined, we assumed a 25% higher concentration than that of Tor1p based on expression data7. The concentration of Tap42p was only available relative to Pph21/22p and Sit4p. We thus indirectly inferred from the above data 1, 4 that the Tip41p concentration exceeds that of Tap42p by about three-fold (Tables S7,S8), and a deviation of ±50% was allowed. Since the degree of phosphorylation of Tap42p in ∆cdc55 and ∆tpd3 mutants was identical 2, the results were considered as one single experiment and an overall concentration of Cdc55p/Tpd3p was used (Table S8).

Table S4 Relative degree of Tap42p phosphorylation in wild type and ∆cd55/∆tpd3 mutants before and after addition of 500 nM rapamycin2. Values were obtained by 32P-labeling of the proteins.

Time after addition of rapamycin [min]
Strain / 0 / 15 / 30 / 60 / 180
Wild type / 100% / 49% / 29% / 12% / 15%
∆cd55/∆tpd3 mutants / 260% / n.d. / n.d. / n.d. / 235%

Table S5 Tap42p-Tip41p complex formation in wild type and ∆sit4 before and after addition of 109 nM rapamycin3. Values were obtained by Adobe Photo Shop 5.5 grey scale analysis of western blots.

Time after addition of rapamycin [min]
Strain / 0 / 30
wild type / 100% / 543%
∆sit4 mutant / 102% / 138%

Table S6 Tap42p-Sit4p complex formation in wild type and ∆tip41 before and after addition of 109 nM rapamycin. Values were obtained by Adobe Photo Shop 5.5 grey scale analysis of western blots.Only the western blot after 30 min was shown, but the values after 60 min are described to be identical 3.

Time after addition of rapamycin [min]
Strain / 0 / 30 / 60
Wild type / 100% / 40% / 40%
∆tip41 mutant / 112% / 82% / 82%

Table S7 The amount of Tap42p, Pph21/22p, Sit4p and Sapp in wild type present in complexes relative to the overall concentration 1, 4.

Complex / Assumed reference values and standard errors 1
Tap42p-Pph21/22p / 10% Tap42p*; M50% / 2% Pph21/22p; M100%
Tap42p-Sit4p / 10% Tap42p; M50% / 5% Sit4p; M100%
Sit4p-Sapp / 60% Sit4p†; M100% / 15% Sapp†; M83%

* 5% Tap42p was found to be present in a complex with Pph21p and 5% in a complex with Pph22p 1

† as stated in4

Table S8 Absolute concentrations of the proteins in the TOR-pathway 5, 6

Protein / Absolute concentration
[# of molecules· cell-1] / Molar concentration **
[nM]
Tor1/2p / 589/723* / 57
Tap42p / ~1’965 / 86
Pph21/22p / 5’620/4’110* / 425
Cdc55/Tpd3p / 8’600/16’900† / 376
Tip41p / 5’710 / 250
Fpr1p / 43’300 / 1’892
Sit4p / 3’970 / 173
Sapp‡ / 11’200 / 489

* Overall concentrations are added because Tor1p/Tor2p and Pph21/Pph22p were considered to have identical functions.

† Both PP2A1 subunits were lumped into a single regulatory subunit, the concentration of which was determined by Cdc55p availability.

‡ As detected for Sap185p

** assuming a cell volume of 3.8•10-11 ml

3.Parameter estimation

For parameter estimation, we employed an unconstrained evolution strategy (ES), more specifically, a custom implementation of the (μ+λ)-ES evolutionary optimization algorithm that considers only bounds on the values of the decision variables8 (kinetic parameters and total protein concentrations). In a comparative study on the identification of kinetic parameters in biochemical networks, the closely related (μ,λ)-ES belonged to the best-performing algorithms9.

For the original model ensemble, since stochastic methods for parameter identification cannot guarantee identification of the global optimum within finite time, we first established appropriate starting points through 100 independent parameter estimations of 440 generations (parents / children ratio of 2:4). From the best initial values, we then searched the parameter space more thoroughly for another 10,000 generations to ensure identification of all those models that are quantitatively consistent with all data.

For the decoupled ensemble (Table S3) and for the ensemble with decoupling and multiple phosphorylation of Tap42p, we employed an iterative estimation algorithm in an outer loop around the individual ES runs (Fig. S1). This iterative approach was performed to enforce uniformity of training errors across each model ensemble and, hence, to increase the reliability of model discrimination because the likelihood of misclassifications due to incomplete parameter estimation would be reduced. Parameter estimation was performed in two ways: (i) either with random seeds as described above or (ii) with a previous solution from a model of the same class (usually the core model) as starting guess. In the latter case, only 2,000 generations were considered for each parameter estimation, in comparison to the random seeds where 10,000 generations followed 100 independent parameter estimations times 440 generations (see above for the original ensemble).


* for random seeds

† Parameter estimation begins here when a previous solution was used as starting guess.

‡ when a previous solution was used as starting guess

Figure S1 Iterative algorithm for parameter estimation used for the decoupled model ensemble. ES: evolutionary strategy; EM: elementary model extensions; θ: set of parameters; εT: mean total model error; σT: standard deviation.

When all of these estimations in iteration k were completed, the total training errors for all elementary models i, i,k(i.k) with the corresponding parameter vector i.k , were used to calculate the errors’ mean and the standard deviation across the model ensemble, εT and σT, respectively. This applied for both starting guesses, i.e. either random seeds or a previous solution, respectively. Solutions for models that fulfilled the error criterion i,k(i.k) ≤ εT ± 1.1 σT, were kept, otherwise another parameter estimation was performed based on the previous run. If the models failed to fulfil the above criterion more than three times, a different initial solution was taken and the parameter estimation for the specific model was started all over again.

Numerical values of estimated kinetic parameters and total protein concentrations are available as Supplementary Data 2 online.The following specifications apply throughout:

-For kinetic parameters, the units are [nM-1·min-1] for association constants and [min-1] for dissociation constants, respectively.

-The upper bound for all parameters during estimation was 103, the lower bound 10-2. With the units for the kinetic parameters this choice of bounds enabled biologically realistic parameter values such as time constants for elementary reactions in the minute to millisecond range. Likewise, these settings enable the estimation of dissociation constants for simple dimerization reactions (KD) between 10-14 M and 10-4 M.

-The association rate constant of rapamycin and Fpr1p (K22) was set between 10-4 and 10-2 to take the diffusion step into account. The dissociation rate constant (K22*) was set equal to 5·K22 (KD = 5 nM, 10).

-Protein concentrations are given in [nM].

-Special modifications of models are indicated by superscripts.

4.Model analysis and validation

4.1Cross-validation

For k-fold cross validation, we chose k = 5 because of the few experimental data available. The data sets were partitioned according to Table S9, taking into account that (i) each partition should contain an approximately equal number of data points, (ii) time-series data was not disrupted, and (iii) each partition should be comparable in terms of the total training error previously estimated and the diversity of data (genotypes, perturbations with rapamycin). Apparently, there exist trade-offs between these objectives, and other partitions could be envisaged.

Table S9 Partitioning of experimental data sets for cross-validation. Detailed descriptions of the experimental data are provided in Tables S4-S8.

Partition / Measured variable / Genotype / +/- Rap. / Number of data points
1 / Tap42p phosphorylation / WT / + / 5
Tip41p phosphorylation / WT / + / 1
2 / Tap42p-Sit4p complex / WT / + / 4
Sit4p-Sapp complex / WT / - / 2
3 / Tap42p-Tip41p complex / WT / + / 3
Tap42p-Sit4p complex / tip41 / + / 3
4 / Tap42p phosphorylation / cdc55 / + / 2
Tap42p-Pph21p complex / WT / - / 2
Tap42p-Tip41p complex / sit4 / + / 3
5 / Tap42p-Pph21p complex / WT / + / 4
Tap42p-Sit4p complex / WT / - / 2

For parameter estimation, we employed the parameter vectors estimated for the respective model ensemble on the entire data set as a starting guess. For each of the k subset of partitions (with k-1 partitions as training data), parameters were estimated for 2,000 generations (2:4 parent / children ratio), and the error in predicting the data in the remaining partition (test set) was determined. The prediction error for each model in the ensemble was determined as the average prediction error for the k test sets.


Figure S2 Prediction errors of the original ensemble in terms of experimental data sets and model classes. Shading denotes relative errors in describing Tap42p~p dephosphorylation (dark grey) and Tap42p-Tip41p complex formation (light grey) in wild type, and for all other experimental data considered (white). All errors are normalized with the maximal prediction error in the two model groups.

Detailed error analysis for the original ensemble revealed that the same experimental data as for the training errors were the main causes for the prediction errors, namely Tap42p~p dephosphorylation and Tap42p-Tip41p complex formation in wild type (Fig. S2). Also, the prediction of Tap42p-Tip41p complex formation constituted the main difference in model performance between elementary extension models 1-7 and the other models in the original ensemble.

Cross-validation provided data to analyze variability in parameter estimation and, thus, assess the suitability of the optimization methods in addition to the independent tests described below (section 4.3). Figure S3 shows average values for the k estimates and the corresponding coefficients of variation for kinetic parameters employed in all original ensemble models. Note that, here, the different data sets used in the individual identifications contribute to variability in the estimated parameter vectors.

Figure S3 Parameter values and variance in parameter estimation. Average parameter values (a, log10 of parameter values) and coefficients of variation (b) for the parameters common to all models in the original ensemble (parameter numbers indicated on the x-axis). For calculating the data, we employed all parameter estimates from five-fold cross-validation.

Fig. S3a shows that for most of the kinetic parameters considered, the estimates varied among the models. It is important to note that, although the same parameters and corresponding reactions are part of all original models, the extensions that lead, among others, to significant differences in model performance, influence these reactions through network effects. Thus, for the core reactions that are highly connected, one can expect a variability of estimates across models. In contrast, the variability of parameter estimates within each model was rather low (Fig. S3b), which points to consistency of the estimation procedure. Moreover, no correlation between the absolute average parameter values and their variability within the models is apparent, which indicates that boundary effects (influence of the constraints on parameter values on the estimates) may have only a minor influence.

4.2Hybrid models


Several of the putative molecular mechanisms described by the 18 elementary model extensions could, in principle, operate simultaneously. To elucidate whether combinations of elementary models would improve description of the experimental data by the original model ensemble, we generated 51 hybrid models that combine the original core model with more than one additional hypothesis specified in Table S1, with a focus on combinations of the seven best elementary extensions (Figure S4). Parameter values for hybrid models were estimated as for the original model ensemble.

Figure S4 Hybrid extensions (x-axis, a) of the original core model obtained by combining elementary models (y-axis, a and b). The pattern of model compositions (a), the corresponding elementary extension models (b) and the normalized training errors for the hybrid extensions (c) and the elementary extensions (d) are shown. The circle in (c) indicates the training error of single-extension model 1.

Only hybrid model 1 showed a marked decrease in the training error compared to the best-performing single-extension model (Fig. S4c). This hybrid model combines kinetic decoupling with multiple phosphorylation of Tap42p. The improved training error for this model that was obtained independently, thus, confirms the results obtained by analysis of the three classes of single-extension models as discussed in the main text. In addition, we could use the parameters of hybrid model 1 (see Supplementary Data 2 online)as suitable starting guesses for the third model ensemble.


In combination with the original model ensemble, we also used the hybrid models to investigate possible correlations between the number of parameter (the degrees of freedom for parameter estimation) and the models’ training errors. A correlation between the two could imply that better performance of a model is mainly due to increased flexibility in describing the few data available and potential overfitting, and not determined by the model structure. We find a weak, statistically significant linear correlation (r2= 0.08, p = 0.02), but Fig. S5 shows that the clustering of the data into approximately two model classes appears to be the dominant feature. This is consistent with the results obtained by cross-validation (relation between prediction error and number of parameters, see main text).

Figure S5 Training errors versus degrees-of-freedom in terms of the number of estimated model parameters for the original model ensemble (filled circles) and the hybrid models (open circles). Training errors were normalized by the maximum error in this set of models. The dashed line shows linear regression.

4.3Analysis of model behaviors

To discriminate models in the original model ensemble according to their behavior, we focused on protein complex formation and protein phosphorylation in the two steady-states before and after addition of rapamycin. The significance of differences in concentration values between model groups Ia, Ib, and II was computed using a two-tailed t-test that probes that two independent samples come from distributions with identical means. Figure S6 summarizes the analysis results for model variables related to Tap42p and Tip41p control.

Figure S6 Comparison of steady-state behaviors between model groups in the original ensemble. States or sums of states relating to Tap42p and Tip41p (complexes indicated by ‘*’, unbound and phosphorylated forms) were considered in steady-state without (‘-‘) or with (‘+’) addition of rapamycin (see grouping on right-hand side). Comparisons that yielded statistically significant differences in model behavior between model groups (assessed by two-tailed t-tests with independent samples, threshold α=0.05) are indicated by grey color. Numbers in the corresponding fields are –log10(p); that is, they represent scores for the significance of the differences.