Model description
We have developed an ordinary differential equation (ODE) model of the MITF-PIAS3-STAT3 system. The graphical representation of the model given in Figure is represented in Systems Biology Graphical Notation as implemented in CellDesigner [
50]. The ODE-model was implemented in MATLAB (Additional file
2 - MATLABcode.zip). In the following the model is presented in terms of the chemical reaction equations and the ordinary differential equations. All protein amounts represented by the state variables in the model are presented in an arbitrary unit (au).
The amount of phosphorylated ERK is represented by a constant (a MAPK signal is represented by an altering of this constant). This kinase is involved in two processes; phosphorylation of RSK1 (1) and phosphorylation of MITF at S73 (2).
In all the equations the subscript p after the protein denotes phosphorylated form. The phosphorylation and de-phosphorylation rate constants are represented by kRp+ and kRp- respectively.
All the equations in (2) represent ERKp catalyzed phosphorylation of MITF at S73. Subscripts after MITF denote phosphorylation state (un-phosphorylated, S73 phosphorylated, S409 phosphorylated and phosphorylated at both sites). Complexes are represented by underscore separated protein abbreviations. Note that the same phosphorylation (kMp73E+) and de-phosphorylation (kMp73-) rate constants are used in all four equations. In addition MITF is auto-phosphorylated at S73 governed by a separate rate constant, kMp73a+:
Further, the phosphorylation of MITF at S409 is catalyzed by the phosphorylated RSK1:
The four different phosphorylation states have different association and dissociation rate constants:
STAT3 is also binding PIAS3:
And finally the phosphorylation of STAT3 is catalyzed by phosphorylated JAK (which is represented by a constant):
A signal propagating through the JAK-STAT3-pathway is represented in the model simulations by an elevation of the constant representing the amount of phosphorylated JAK.
These chemical reactions were translated into an ODE model by anticipating simple mass action reaction kinetics. See Table for a comprehensive list of parameters. Equation (11), describing the dynamics of unbound and un-phosphorylated MITF, is provided in detail below. One differential equation for each phosphorylation and complex state of each protein is given. These equations are in the same way determined from the chemical equations and are therefore only commented when features differ.
In this model, MITF can be present in 8 different states and is therefore represented by 8 differential equations, each containing parts reflecting the chemical reactions involving that particular state.
In differential equation (11), where the rate of change of un-phosphorylated and un-bound MITF is described, PMITF is the production rate, kMp73E+[MITF][ERKp] is the rate of MITF leaving this state because it gets phosphorylated at S73 (reaction equation (2)), kMp73-[MITFp73] is the rate of MITF entering this state from the S73 single phosphorylated state due to de-phosphorylation (reaction equation (2)), kMp73a+[MITF] is the rate of MITF leaving this state because of auto-phosphorylation at S73 (reaction equation (3)), kMp409+[MITF][RSK1p] is the rate of MITF leaving this state because of phosphorylation at S409 (reaction equation (4)), kMp409-[MITFp409] is the rate of MITF entering this state from the S409 single phosphorylated MITF state because of de-phosphorylation (reaction equation (4)), kMass[MITF][PIAS3] is the rate of MITF leaving this state because of MITF-PIAS3 association (reaction equation (5)) and kMdiss[PIAS3_MITF] is the rate of MITF entering this state due to MITF-PIAS3 dissociation (reaction equation (5)). The last term γMITF[MITF], represents the MITF degradation. The total amount of RSK1 is kept constant, while its phosphorylation state is dynamically determined by equations (12) and (13).
The dynamics of S73 phosphorylated MITF (MITFp73) is given in equation (14)
Where γMITFp73 is the degradation rate that applies to the ubiquitinated ratio R of the amount of MITFp73. In the dynamics of S409 phosphorylated MITF,
simple linear degradation is applied. When MITF is phosphorylated on both sites, the sum of the two degradation rates applies:
Further, these four phosphorylation states of MITF bind PIAS3, each with their own association and dissociation rate constants, yielding differential equation (17) for PIAS3, (18) for the MITF-PIAS3-complex, (19) for the S73 phosphorylated MITF-PIAS3-complex, (20) for the S409 phosphorylated MITF-PIAS3-complex and (21) for the double phosphorylated MITF-PIAS3-complex.
PIAS3 has a linear degradation.
STAT3 is binding PIAS3 while phosphorylated. The dynamics of the PIAS3-STAT3-complex is given in equation (22) and un-phosphorylated and phosphorylated STAT3 in equation (23) and (24) respectively.
To represent higher stability by proteins in complexes, the degradation rate of the complex is set to 20% of the mean of the degradation rates assigned to the constituting proteins in un-bound form.
S73 phosphorylated MITF get ubiquitin tagged for degradation by phosphorylated ERK. Introducing another modification site would cause another doubling of the number of state variables representing MITF. To reduce complexity, we chose to represent the ubiquitination with one differential equation (25), describing the ratio R of the total amount of S73 phosphorylated MITF that also is ubiquitinated.
Where A is defined as
The three parts of equation (25) correspond to the three processes altering the ratio R: (i) The ubiquitination of un-ubiquitinated MITF, (ii) the ubiquitination dependent degradation of ubiquitinated MITF and (iii) the production of new un-ubiquitinated MITF. In (i), the ubiquitination rate constant ku, denotes the rate of change of R and thus have the unit au-1 min-1. To derive part (ii), let R' be the ratio R after one minute of ubiquitin mediated degradation. Then the rate of change due to this effect is -(R - R') per minute, where
and B = [MITFp73] + [MITFpp]. The term a = 1 and has the unit minute and does not affect the calculation, but is needed for the correct unit representation. Part (iii) is derived in a similar way: Let R'' be the ratio R after one minute of MITF production. Then the rate of change due to this effect is -(R - R'') per minute, where
MITF activity
The readouts of the model are the amounts of the different phosphorylation and complex states of the proteins involved. To be able to view this model in the light of available experimental results, these levels have to be interpreted in a way comparable to the lab-generated data. The transcriptional activity of MITF and STAT3 is measured either as luciferase activity of a transfected reporter gene or as mRNA-levels of downstream genes.
We thus devised a function f([MITF],[MITFp73],[MITFp409],[MITFpp]) to map the amounts of the different phosphorylation states of MITF to a downstream transcriptional activity. The production of any gene product can be modelled by a function f(S1,S2,....,SN), where Sj is the concentration of the N species affecting the transcription and translation of the particular gene product. This function can be described at various levels of precision from the basic principles of chemical kinetics, via Michaelis-Menten enzymatic kinetics to low order polynomials that function as an ansatz for the actual molecular mechanisms. The function that maps from the amounts of the MITF phosphorylation states to production rate of downstream genes will probably differ between target genes on a higher resolution level. However, for the resolution level of the current effort, one mapping function is sufficient. We have approximated the MITF activity by a sum of first order polynomials, which in effect is a weighted sum:
Values for the transcriptional activity of different MITF mutants are given in [
14]. From these we have calculated the factors to be
AM = -0.11,
AM73 = 0.44,
AM409 = 0.11 and
AMpp = 0.56. This simple representation of MITF transcriptional activity is able to take negative values, which does not make sense biologically; therefore,
A0 is set to the relative high arbitrary value 10, which is sufficiently high to avoid negative MITF activity in this study.
Parameterizing the model
While this model contains a large number of parameters, we can make at least a reasonable order of magnitude estimate for most of them from empirical biological data. In the experiments used for model development and parameter fitting, neither implicit nor explicit absolute protein amounts have been provided. The cells are either transfected, activated or both, and thus the input protein amounts and perturbations will by this fact be qualitative. The results are presented as qualitative protein or complex amounts (Western blot) or transcriptional activity for MITF and PIAS3 (measured by luciferase assays or PCR-measurements of target genes mRNA levels). These results are quantified relative to a control, but not in absolute terms.
An estimate for the degradation rates of PIAS3 and STAT3 is calculated directly from the protein stability index (PSI) as presented in [
25] (
γSTAT3 = 0.002 and
γPIAS3 = 0.008). The PSI derived degradation rate for MITF is 0.008, but MITF is known to be relatively stable in un-phosphorylated form while ERK mediated phosphorylation on S73 also tags MITF for proteasome mediated degradation [
13,
14]. To match degradation time series data [
13,
14] MITF was assigned a low degradation rate in the un-phosphorylated state,
γMITF = 0.0012, and relatively high degradation rate when phosphorylated on S409 (regardless of phosphorylation status on S73),
γMITFp73 = 0.01. The proportion of MITF phosphorylated on S73 (regardless of phosphorylation status on S409) is degraded only when ubiquitinated (and then with a rather high degradation rate,
γMITFp73 = 0.02). The ubiquitination rate constant is adjusted to
ku = 0.0001 to yield degradation in accordance with data [
13,
14]. When the proteins are in complex, they are thought to be more stable than as free molecules. As an approximation the degradation rates of all complexes were set to 20% of the average of the degradation rate of the two components.
The production rates for MITF, PIAS3 and STAT3 are complex functions of the concentration of a large number of substances, of which none are involved in this model. Therefore these production rates are represented here by constants whose values are relative to the corresponding degradation rates determined by the steady state levels of each of the three proteins. It is important for the proper function of this module that these three proteins are present at approximately equimolar levels. The values are therefore set to
pMITF = 1,
pPLAS3 = 1.262 and
pSTAT3 = 0.211 to meet this criterion. This yields a steady state amount of approximately 100 (au). The level of phosphorylated ERK in the resting cell is set to
ERKp = 10 and is elevated to approximately 1000 whenever a MAPK signal is simulated. The total amount of RSK1 is set to 500 of which the phosphorylation is determined dynamically (equation (12) and (13)). The phosphorylation speed for RSK1 is found in [
14], and the phosphorylation and de-phosphorylation rate constants calculated to
kRp+ = 0.0004 and
kRp- = 0.04, respectively. The Western blot analysis of MITF presented in [
13] and [
14] has one band interpreted as un-phosphorylated MITF and one band interpreted as phosphorylated MITF. No interpretation is provided on which phosphorylation states this latter band may represent. By interpreting this band as the sum of all phosphorylated MITF, we have adjusted the phosphorylation rate constants of MITF to make the model behaviour fit the experimental data. The rate constants are set to
kMp73E+ = 0.00015,
kMp73- = 0.03,
kMp73a+ = 0.025,
kMp409+ = 0.0001 and
kMp409- = 0.04. With no direct measurements of the STAT3 phosphorylation rate, we assigned the same values to the STAT3 phosphorylation/de-phosphorylation rate constants as we did for RSK1:
kSp+ = 0.0002 and
kSp- = 0.04. The constant representing phosphorylated JAK is set to 10 in the un-activated state and approximately 1000 in the activated state.
With no direct measurements of the MITF-PIAS3 and STAT3-PIAS3 association and dissociation rate constants, they were inferred from qualitative results and assumptions. We anticipated that this process is faster than the phosphorylation process, which means that both association and dissociation rate constants should be relatively high. MITF has lower affinity to PIAS3 when phosphorylated at S409 and higher affinity when phosphorylated at S73, compared to un-phosphorylated MITF [
20]. The double phosphorylated MITF exhibits an intermediate affinity comparable to that of un-phosphorylated MITF [
20]. When the MITF-PIAS3 association/dissociation rate constants were adjusted so that the model emulated the results from [
20], they were assigned the following values:
kMiss = 0.01,
kMdiss = 1,
kMp73ass = 0.03,
kMp73diss = 0.5,
kMp409ass = 0.0001,
kMp409diss = 1,
kMppass = 0.01 and
kMppdss = 1. The STAT3-PIAS3 association/dissociation rate constants were adjusted in the same way to data from [
22] and are assigned the values
kSpass = 0.005 and
kSpass = 0.2.
Formalizing the experiments
We extracted a set of representative experiments from four relevant publications [
14,
19,
20,
22] (Table ). The numbering in the table provided a unique ID for each experiment, and was used throughout the article. The perturbations of the cells that were reported in the published experiments, like transfection of mutated genes or receptor activation, were assigned digital counterparts in the model through alterations of the input values, starting values, or some of the parameters. In the following, these details, as well as the criteria used for the sensitivity analysis, are given for each of the experiments.
Experiment # 1 is a simulation of the experiment presented in [
14] Figure . Here, the authors use kinase assays to investigate the temporal development of ERK and RSK1 kinase activity after stimulation of melanoma cells. The stimulation was simulated by elevation of the amount of phosphorylated ERK to 1000 from the default value of 10. The success criterion used in the sensitivity analysis was that the phosphorylated RSK1 level should reach 90% of its maximum between 3 and 10 minutes after activation.
Experiment #2 simulates the pre-stimulation distribution among the MITF phosphorylation states in the experiments presented in [
14], Figure and . To emulate the growing cell culture the level of phosphorylated ERK was elevated from the default value of 10 to the intermediate level of 80. Steady state values were found for the state variables and the total amount of phosphorylated and un-phosphorylated MITF was compared. The success criterion used in the sensitivity analysis was that both categories should be between 25% and 75% of the total amount.
Experiment #3 is also based on [
14], Figure and . Here, the phosphorylation state of MITF after 30 minutes of activation is considered. The cell stimulation is simulated by elevation of the levels of phosphorylated ERK and JAK from the default value of 10 to 1000. The success criterion used in the sensitivity analysis was that more than 80% of the MITF is phosphorylated after 30 minutes.
Experiment #4 is also based on [
14], Figure and . Here, the degradation profile is considered. The cell stimulation is simulated by elevation of the levels of phosphorylated ERK and JAK from the default value of 10 to 950. The success criterion used in the sensitivity analysis is that less than 50% of the MITF is degraded after 1 h and at least two thirds are degraded after 5 h.
Experiment #5 is a simulation of the experiment presented in [
19], Figure . Here, the authors have transfected MITF, a MITF reporter gene, and different amounts of PIAS3 into NIH 3T3 fibroblasts (endogenously not expressing MITF). To mimic MITF transfection, the MITF production rate is increased from 1 to 18. To mimic the four different amounts of PIAS3 transfection, the PIAS3 production rate was increased from 0.211 to 0.5, 1, 2 and 4. The model was simulated for 2880 minutes (2 days) for each perturbation and the resulting MITF activity was compared for the four PIAS3 levels. The success criterion used in the sensitivity analysis was that MITF activity should be monotonically decreasing with increasing PIAS3 and that the MITF activity at the highest PIAS3 amount should be less than the half that of the lowest PIAS3 amount.
In experiments #6 and #7, we have mimicked the experiments presented in [
20], Figure . Here the authors have investigated the association between MITF and PIAS3 in BL6-B16 cells before and after 10 and 30 minutes of activation with c-kit ligand. In experiment #6 only GFP-PIAS3 was transfected, while in experiment #7, the cells were transfected with both MITF and PIAS3. In experiment #6, activation is performed by increase of the levels of phosphorylated ERK and JAK from 10 to 500 and 50, respectively. In experiment #7, transfection was simulated by elevation of the MITF production rate from 1 to 10 and the PIAS3 production rate from 1.262 to 10 and simulated for 2880 minutes (2 days). Thereafter, activation was simulated by increase of the amount of phosphorylated ERK and JAK from 10 to 1200 and 250, respectively. The success criterion used for the sensitivity analysis in both experiments was that the total amount of MITF-PIAS3-complex should be higher after 10 minutes of stimulation compared to before stimulation and that the amount of complex should be lower after 30 minutes of stimulation than after 10 minutes.
Experiment #8 simulated the experiment from [
20], Figure , left. The authors co-transfected NIH 3T3 cells using MITF, with a constitutively active RSK1 plasmid and with a reporter gene to read MITF activity. To simulate transfection of MITF, the MITF production rate was elevated from 1 to 7. The amount of RSK1 is not dynamically determined in the model, and thus the transfection of this product is approximated with an elevation of the static amount from 500 to 5000. The manipulated model was simulated for 2880 minutes (2 days), and the MITF activities with and without RSK1 were compared. The success criteria used in the sensitivity analysis was that the two cases should not differ more than twofold.
Experiments #9 to #12 represent simulations of the experiment presented in [
20] Figure . Here, the authors transfected NIH 3T3 cells with MITF or the MITF-S409A mutant, with constitutively active RSK1 and with a reporter gene to monitor MITF transcriptional activity. To simulate the transfection of MITF, the MITF production rate was elevated from 1 to 7. To simulate the transfection of PIAS3, the PIAS3 production rate was elevated from 1.262 to 4. The amount of RSK1 is not dynamically determined in the model, and thus the transfection of this product is approximated with an elevation of the static amount from 500 to 5000. The activation level was elevated from 10 to 20 to emulate the activity in growing cells. After simulation for 2880 minutes (2 days), the MITF activities, with and without PIAS3 were compared for each of the four experiments. For each experiment the manipulations were as follows: without RSK1, with wild-type MITF (experiment #9), with RSK1 and wild-type MITF (experiment #10), without RSK1, with mutant MITF (experiment #10) and with RSK1, with mutant MITF (experiment #11). The success criterion used in the sensitivity analysis was for each experiment: Experiment #9: MITF activity in the case with PIAS3 is less than 42% and greater than 2% of that without PIAS3. Experiment #10: MITF activity in the case of PIAS3 is less than 84% and greater than 44% of that without PIAS3. Experiment #11: MITF activity in the case of PIAS3 is less than 55% and greater than 15% of that without PIAS3. Experiment #12: MITF activity in the case of PIAS3 is less than 50% and greater than 10% of that without PIAS3.
Experiment #13 to #16 are simulations of the experiment presented in [
20] Figure . In this experiment the authors investigate the inhibition of transcriptional activity of S/D mutants of MITF by PIAS3. This is achieved by transfection of wild-type MITF (experiment #13), MITF-S73D (experiment #14), MITF-S409D (experiment # 15) or MITF-S73/409D (experiment #16) and PIAS3 and a reporter gene into NIH 3T3 cells. We simulated the MITF transfection by elevation of the MITF production rate from 1 to 1.5, and the PIAS3 transfection was simulated by elevation of the PIAS3 production rate from 1.262 to 4.5. Here we have also set all MITF levels to zero before starting the simulation, since NIH 3T3 fibroblasts do not express endogenous MITF. To simulate MITF-S73D mutation we have set the MITF S73 de-phosphorylation rate constant to zero while the MITF S73 auto-phosphorylation rate constant is elevated from 0.025 to 5. To simulate MITF-S409D mutation, the RSK1 phosphorylation rate constant is elevated form 0.0004 to 0.04, the RSK1 de-phosphorylation rate constant is decreased from 0.04 to 0.004, the MITF S409 phosphorylation rate constant is elevated from 0.0001 to 0.01 and the MITF S409 de-phosphorylation rate constant is set to zero. The model was simulated for 2880 minutes (2 days) and the MITF transcriptional activity with PIAS3 was compared with the case without PIAS3 for wild-type MITF, and for each mutation form. The success criterion used for the sensitivity analysis was for each experiment: Experiment #13: MITF activity in the case of PIAS3 is less than 70% and greater than 30% of that without PIAS3. Experiment #14: MITF activity in the case of PIAS3 is less than 50% and greater than 10% of that without PIAS3. Experiment #15: MITF activity in the case of PIAS3 is less than 115% and greater than 75% of that without PIAS3. Experiment #16: MITF activity in the case of PIAS3 is less than 85% and greater than 45% of that without PIAS3.
Experiments #17 and #18 are simulations of the experiment presented in [
22] Figure . Here, the authors investigate PIAS3-STAT3 association in response to stimulation of NIH 3T3 cells transfected with STAT3, PIAS3 and MITF or the MITF S409A-mutant. The transfection is simulated by elevation of the MITF production rate from 1 to 10, the PIAS3 production rate from 1.262 to 10 and the STAT3 production rate from 0.211 to 10. The MITF S409A mutation was simulated by setting the MITF S409 phosphorylation rate constant to zero. The activation was simulated by elevation of the amount of phosphorylated ERK and JAK from 10 to 1000. In experiment #17, wild-type MITF, PIAS3 and STAT3 were all transfected and simulated for 2880 minutes (2 days), thereafter simulated for 15 more minutes with and without stimulation. The amount of PIAS3-STAT3-complex was compared. The success criterion used in the sensitivity analysis was that the amount of PIAS3-STAT3-complex should be at least 25% higher with stimulation compared to the case without stimulation. In experiment #18, S409A MITF, PIAS3 and STAT3 were transfected and simulated for 2880 minutes (2 days), and thereafter simulated for 15 more minutes with or without stimulation. The amount of PIAS3-STAT3-complex was compared with the stimulated case in experiment #17. The success criterion used in the sensitivity analysis was that the amount of PIAS3-STAT3-complex in the stimulated case in experiment #17 should be at least 25% higher than the amount of PIAS3-STAT3-complex in any of the two cases in #18.
Experiments #19 to #24 are simulations of the experiment presented in [
22] Figure . Here, the authors have investigated STAT3 transcriptional activity as a response to activation, transfection of PIAS3, and transfection of various amounts of wild-type and S409A mutated MITF in NIH T3T cells. The activation was simulated by elevation of the amounts of phosphorylated ERK and JAK from 10 to 1000. The transfections of PIAS3, STAT3 and the two different MITF amounts were simulated by elevation of the PIAS3 production rate from 1.262 to 7, the STAT3 production rate from 0.262 to 5, and the MITF production rate from 1 to 10 or 50, respectively. The S409A mutation was simulated by setting the MITF phosphorylation rate constant to zero. The model was run for 2880 minutes to simulate the incubation and then for 360 minutes to simulate the activation. The background luciferase activity of STAT3 was determined by elevation of the STAT3 production rate and running a model simulation for 3240 (2880+360) minutes. The amount of phosphorylated STAT3 in the end was interpreted as an estimate for the luciferase activity of STAT3 without activation and without PIAS3 or MITF. For experiments #19 to #24, the amount of phosphorylated STAT3 was compared to this background level. In experiment #19, the cells were transfected with STAT3 and activated. The success criterion used in the sensitivity analysis was the level of phosphorylated STAT3 between 10 and 20 times the background level. In experiment #20, the cells were transfected with STAT3 and PIAS3 and activated. The success criterion used in the sensitivity analysis was the level of phosphorylated STAT3 between 2.67 and 5.33 times the background level. In experiment #21 the cells were transfected with STAT3, PIAS3 and the smaller amount of MITF and activated. The success criterion used in the sensitivity analysis was the level of phosphorylated STAT3 between 2.67 and 5.33 times the background level. In experiment #22 the cells were transfected with STAT3, PIAS3 and the larger amount of MITF and activated. The success criterion used in the sensitivity analysis was the level of phosphorylated STAT3 between 6.67 and 13.33 times the background level. In experiment #23 the cells were transfected with STAT3, PIAS3 and the smaller amount of S409A mutated MITF and activated. The success criterion used in the sensitivity analysis was the level of phosphorylated STAT3 between 5.67 and 11.33 times the background level. In experiment #24 the cells were transfected with STAT3, PIAS3 and the larger amount of S409A mutated MITF. The success criterion used in the sensitivity analysis was the level of phosphorylated STAT3 between 10 and 20 times the background level.
Experiment #25 simulated the experiment presented in [
22] Figure . Here, the authors have investigated MITF transcriptional activity in response to transfection of PIAS3, STAT3-Y705F mutant and a constitutively active STAT3 mutant (STAT3-C) in NIH 3T3 cells. The transfections were simulated by elevation of the production rate of MITF from 1 to 2, PIAS3 from 1.262 to 2 and STAT3 from 0.211 to 20. The NIH 3T3 cells were simulated by decreasing the starting values of all MITF and PIAS3 states to 1 and in simulations where MITF or PIAS3 were not transfected; their production rate were set to zero. The STAT3-Y705F mutant was simulated by setting the STAT3 phosphorylation rate constant to zero. STAT3-C was simulated by increase the STAT3 phosphorylation rate constant from 0.0002 to 0.02 and setting the de-phosphorylation rate constant to zero. Six different simulations were performed and the MITF transcriptional activity was compared: (i) Cells transfected with STAT3-Y705F, (ii) cells transfected with STAT3-C, (iii) cells transfected with MITF and STAT3-Y705F, (iv) cells transfected with MITF and STAT3-C, (v) cells transfected with MITF, PIAS3 and STAT3-Y705F and (vi) cells transfected with MITF, PIAS3 and STAT3-Y705F. The success criterion used in the sensitivity analysis was that the MITF transcriptional activity from simulation (i) should be between 50% and 150% of the MITF activity form (ii) and the MITF activity from (iii) should be between 80% and 120% of the MITF activity from (iv) and the MITF activity from (v) should be between 12.5% and 100% of the MITF activity from (vi) and the MITF activity from (i) should be less than the MITF activity from (iii) and the MITF activity from (iii) should be greater than 200% of the MITF activity from (v).
Experiment #26 simulates the experiment presented in [
22] Figure and . Here, the authors use mRNA levels of MITF and STAT3 target genes to investigate MITF and STAT3 transcriptional activity after stimulation of mast cells derived from wild type or mutant mice. The MITF protein from the mutated mice lacks its ability to bind PIAS3. The activation is simulated by elevation of the amounts of phosphorylated ERK and JAK from 10 to 1000. The mutation is simulated by setting all four MITF-PIAS3 association rate constants to zero. The success criterion used in the sensitivity analysis was for cells from wild type mice: MITF activity after 30 minutes should be more than twice the MITF activity without stimulation, STAT3 activity at 4 hours should be more than twice the STAT3 activity without stimulation and for cells from mutated mice: MITF activity at 30 minutes and at 4 hours should be less than twice the MITF activity without stimulation and STAT3 activity at 30 minutes and 4 hours should be less than twice the STAT3 activity without stimulation.
Experiments #27 and #28 are simulations of the experiments presented in [
22] Figure and . In experiment #27 the authors have investigated the STAT3 transcriptional activity in response to transfection of wild type MITF and a MITF mutant emulating MITF phosphorylated at S409. In experiment #28 the authors have investigated MITF transcriptional activity in response to transfection of constitutively active STAT3 (STAT3-C). The transfections were simulated by elevation of the STAT3 production rate from 0.211 to 5 and the MITF production rate from 1 to 5. The MITF mutant was simulated by setting the MITF-S409 phosphorylation rate constant to 5 and the MITF-S409 de-phosphorylation rate constant to zero. The STAT3-C mutant was simulated by increasing the STAT3 phosphorylation rate by 25 times and setting the de-phosphorylation rate to zero. The model was simulated for 2880 minutes without and with each of the three transfections. The success criterion used in the sensitivity analysis was for experiment #27 that the STAT3 activity had increased when the cells were transfected with MITF or the MITF mutant compared to without transfection and that the increase was higher when cells were transfected with wild type MITF compared to mutated MITF. The success criterion used in the sensitivity analysis for experiment #28 was MITF activity increased by at least 10%.
Sensitivity analysis
The ability of the model to mimic the experiments described above was tested for values of the core parameters in the multidimensional neighbourhood of the default values. Each parameter was sampled from a distribution being uniform on the logarithmic scale, with a sampling area reaching from one half to double the default value. For each set of values of the core parameters, the experiments were simulated and the result (success or failure) was recorded for each experiment. This procedure was repeated 106 times resulting in 106 parameter sets, each with a corresponding result vector recording success or failure for each experiment.
In order to study the sensitivity of each experiment
i's success rate to variation in a parameter
j, we sorted all parameter sets according to parameter
j and divided them into 100 bins of 10000 parameter sets each. We used
fi , the o verall success rate of experiment
i (see Table , column 6) and

, the success rate in bin
k to computed the sensitivity measure

. In order to remove clearly insensitive parameters, we set to zero all
si, j that did not exceed the maximum value observed under 10000 permutations of the parameter values.