Model calibration

6.5.1 Experimental data for model calibration

Various stimulation experiments were performed and several pathway components were monitored by Peter Nickel in primary hepatocytes as a premise for data-based modelling.

The following five stimulation conditions were tested experimentally: (i) standard TGFβ stimulation (1 ng/ml), (ii) weak TGFβ stimulation (0.05 ng/ml), (iii) standard TGFβ stimulation in the presence of the generic transcriptional inhibitor actinomycin D, (iv) standard TGFβ stimulation with addition of the specific TGFβ receptor kinase inhibitor SB431542 after 50 min, and (v) standard TGFβ stimulation with the ligand being removed after 50 min by medium exchange. Being central to TGFβ signalling, total SnoN, total Smad2, Smad2 phosphorylation and the amount of Smad4 co-immunoprecipitated with Smad2/3 were measured in whole cell lysates. Total and phosphorylated Smad2 pools were additionally measured in cytosolic and nuclear extracts to get insights into the dynamics of nucleocytoplasmic shuttling. Moreover, the amount of processed TGFβ in the cell culture medium was measured in the presence and absence of actinomycin D. All time course data sets used to calibrate the kinetic model by global parameter estimation are summarised in Table 6.1. For background correction, the lowest value of a time course was generally subtracted from all data points for all observed species except for total Smad2, for total SnoN and for TGFβ measured in the cell culture medium.

Fig. 6.7 Model fits to nuclear and cytoplasmic Smad2/4 data.

The kinetic parameters of the model were estimated by fitting the model to experimental data (data points). Solid lines represent the corresponding simulated time series obtained for the best-fit parameters. Relative concentrations of phosphorylated and total Smad2 (A) in cytoplasmic extracts after stimulation with 1 ng/ml TGFβ in absence (black) or presence (grey) of actinomycin D and (B) in cytoplasmic and nuclear extracts after stimulation of cells with 1 ng/ml TGFβ were determined by quantitative immunoblotting. Spline-based error estimates were used (Section 6.5.2).

Cytoplasmic Smad signalling +/- Actinomycin D

Cytoplasmic Smad signalling (standard stimulation)

Nuclear Smad signalling (standard stimulation)

A

B

6.5.2 Scaling factors and error estimation

Scaling Factors: Western blot measurements are typically semi-quantitative in the sense that the results can only be expressed in arbitrary units, but not in absolute concentrations. Thus, scaling factors were introduced to fit the model (formulated in absolute concentrations) to time course data given in arbitrary units. An individual scaling factor was included for each observable (total Smad2, phosphorylated Smad2, Smad4 (CoIP), total SnoN, SnoN (CoIP), TGFβ) in the respective compartments, as indicated in Table 6.1, where each box corresponds to a single scaling factor. It is important to note that a single scaling factor was used for different stimulation conditions and for SnoN-depleted vs. wildtype cells in order to take biologically relevant scaling into account (Table 6.1).

Fig. 6.8 Model fits to whole-cellular SnoN data under various stimulation conditions.

The kinetic parameters of the model were estimated by fitting the model to experimental data (data points). Solid lines represent the corresponding simulated time series obtained for the best-fit parameters. Cells were either stimulated with (A) 1 ng/ml TGFβ or (B) with 0.05 ng/ml TGFβ or (C) with 1 ng/ml TGFβ and actinomycin D was added at 0 min or (D) SB431542 was added after 50 min. Relative concentrations of total SnoN in whole cell lysates were determined by quantitative immunoblotting and are directly comparable between conditions (A-D) due to common scaling. Error estimation was based on a spline-derived error model (Section 6.5.2).

The quality of a mathematical model is determined by its deviation from the data relative to the measurement error, so that the variability of the experimental data had to be estimated.

All measurements listed were performed at least twice with similar results. As replica measurements frequently differed in the observed time points, errors could not be directly derived by calculating standard deviations. The errors of 14 densely measured time courses (each comprising 18 or more time points) were therefore estimated by taking the distance of each data point from a cubic spline. From these spline-based error estimates linear error models were derived for the observables total Smad2, phosphorylated Smad2, Smad4 (CoIP), SnoN (CoIP) and TGFβ. For each time series y(ti), these error models, considering relative errors E_rel and absolute errors E_abs, had the form

E_tot,i = E_rel ⋅ |y(ti)| + E_abs ⋅ C.

Here, the constant C is given by the 0.9-Quantile of the whole time course, i.e., C = Q.9(y(t1...ti)). These error models were used to estimate the errors in less densely sampled time courses (Table 6.1). For total SnoN, no spline-based error model could be derived from the measurements, so that typical values (Erel = 10%; Eabs = 15%) were applied. The relevance of SnoN-mediated feedback was confirmed by analyzing Smad signalling in cells harbouring reduced levels of functional SnoN (see below). The corresponding experiments in SnoN siRNA cells and in knock-in cells expressing mutant SnoN were performed in triplicates with identical time points analysed in each experiment, so that the errors could be estimated by calculating standard deviations. For a comprehensive overview of the error estimation for different data sets see Table 6.1.

6.5.3 Parameter estimation

Mathematical modelling was done using the multi-experiment fitting Matlab toolbox PottersWheel (www.potterswheel.de, Maiwald et al., under revision). The reaction scheme (Fig. 6.4) was translated into a set of 19 ordinary differential equations. The final model additionally included 9 auxiliary variables (see Equations F.1 and F.2), which were required to model a delay in SnoN induction and the mRNA decay in response to actinomycin D treatment. Several of the 57 model reactions were assumed to proceed with identical kinetics, so that 36 kinetic parameters were subject to fitting. All parameters were optimised globally, that is, a single parameter set was able to describe all time courses. During the parameter estimation, each parameter was restricted to the range specified in Tables F.1 – F.3. According to Table 6.1, there are 10 scaling parameters, all of which were also fitted, and allowed to vary in a large range (10^-5 – 10⁵). The initial conditions (i.e., the initial protein concentrations) are not subject to separate fitting, as they are determined by the kinetic parameters for basally expressed species (see Section F.3) or, in the case of the extracellular TGFβ concentration, by the experimental conditions (Section 6.5.4). The start values for the auxiliary mRNA species were also fixed and set to 1 (see Appendix F). Finally, the experimentally measured volumes of medium, cytosol and nucleus were set to experimentally measured values (see Section 6.4.1), and were also excluded from the fitting procedure. PottersWheel uses the weighted χ² value for parameter optimisation,

N 2

2 Model Meas

2

i=1 Meas

(y (i) - y (i))

= (i)

χ ∑ σ

N is the number of all data points and σ_Meas(i) is the error of the i-th measurement, y_Meas(i). In order to circumvent local minima, a two-step strategy was applied. First, 2800 quasi randomly distributed positions in the space of physiologically reasonable parameter values (see Supplementary Tables F.1 – F.3) were used as starting conditions. Fitting was applied with a deterministic trust region optimiser with a χ² tolerance of 10^-4, a fit parameter tolerance of 10^-4, and a maximum of 600 iterations. Then, the best fit was disturbed and refitted 1600 times, i.e., the fitted parameter values were transformed to p_new = p_old × 10^{(s × e)}, with s = 0.3 and e being normally distributed with variance 1 and mean 0.

6.5.4 Input functions

The model trajectories depend on the characteristics of five external input functions, termed u1 – u5 in the differential equations (Equations F.1 and F.2). The input u1 equals the initial TGFβ concentration, and thus determines the strength of stimulation. Ligand removal experiments by medium exchange were simulated by changing the input u₅ from 1 to 0, which abolishes de novo ligand binding to TGFβ receptors in the model (see v7 in Equations F.2). The input u₃ controls receptor-mediated Smad2 phosphorylation (v₁₆ in Equations F.2),

and was changed from 0 to 1 to simulate addition of the TGFβ receptor kinase inhibitor SB431542. In some experiments, TGFβ signalling was monitored over 10h in the presence of the general transcription inhibitor, actinomycin D, and this was modelled by setting the input u₂ = 1 (with u₂ = 0 otherwise). The setting u₂ = 1 induces an exponential decay of the auxiliary mRNA species (v₁, v₃, v₅, v₄₃ in Equations F.2) and inhibits Smad-induced SnoN synthesis (v₄₅ in Equations F.2), thereby mimicking general inhibition of transcription in the model. Finally, experiments in SnoN-depleted hepatocytes were simulated by setting u₄ = 1 (with u₄ = 0 otherwise), as this blocks both induced and basal SnoN expression (v₄₄ and v₄₅ in Equations F.2). The parameter k₁ controls the effective degree of SnoN depletion, and was set to 1 in experiments with the Smad binding-deficient SnoN mutant (complete depletion), while the value was adjusted to 0.7 to simulate siRNA-mediated SnoN knock-down (partial depletion).

6.5.5 Results of model calibration

As shown in Figs. 6.5 – 6.8, the model could simultaneously be fitted convincingly to all data sets. The final best fit had a χ² value of 355. Given N = 506 data points, the model hypothesis could not be rejected (p < 0.05).

The amount of active TGFβ in the medium continuously decreases over the 10 h time course (Fig. 6.5A). This seems to be due to receptor-mediated ligand endocytosis and degradation, as the observed ligand decay can be quantitatively described by the model if experimentally measured kinetic parameters are assumed for the receptor trafficking module (Table F.1).

Interestingly, the measured ligand decay time course is not affected by actinomycin D treatment (Fig. 6.5A), thus suggesting autocrine stimulation does not play a significant role in primary mouse hepatocytes, even though induction of TGFβ and of ligand processing enzymes was observed at the mRNA level (Fig. 6.2).

Despite the strong decay in extracellular TGFβ levels, sustained Smad2 phosphorylation and Smad2-Smad4 complex formation occurred in primary mouse hepatocytes (Fig. 6.5B and C), and the model suggests that is due to strong saturation at the receptor level. Incubation of cells with actinomycin D induced some degradation of Smad2 (Fig. 6.5D). Yet, the level of phosphorylated Smad2 was not significantly affected by actinomycin D treatment (Fig. 6.5B), while the amount of Smad4 co-immunoprecipitated to Smad2 dropped to approximately half of the control (Fig. 6.5C). Thus, Smad2 phosphorylation and Smad2-Smad4 complex show distinct responses towards actinomycin D, while they behave very similarly under various other stimulation protocols (Fig. 6.6). This suggests that complex formation between Smad2 and Smad4 (but not Smad2 phosphorylation) might be controlled by transcriptional feedback loops, as analysed in more detail below (Section 6.6).

Immunoblot measurements of Smad signalling in nuclear and cytoplasmic compartments, respectively, revealed weak translocation of Smad proteins into the nucleus in TGFβ-stimulated primary mouse hepatocytes (Fig. 6.6). In the model, the relatively small stimulus-induced change in both, nuclear and cytoplasmic, Smad2 pools could only be explained if Smad2 was roughly equally distributed between both compartments in resting cells. In line, it has recently been reported in the literature that Smad2 is present at similar amounts in the cytoplasm and in the nucleus of unstimulated primary mouse hepatocytes [364] .

The SnoN protein expression level was measured under various stimulation conditions (Fig.

6.8). SnoN was induced with a delay of ~ 60 min in primary mouse hepatocytes (Fig. 6.8A), as reported previously for cell lines [343] . Actinomycin D incubation led to a rapid decline in SnoN protein levels (Fig. 6.8C). In the model this was mainly due to TGFβ-induced SnoN degradation mediated by phosphorylated Smad2 (see Section 6.4.4) which cannot be

compensated for by TGFβ-induced SnoN synthesis in actinomycin- treated cells. The alternative

Fig. 6.9 Experimental testing of model predictions for mutant SnoN form abrogating SnoN-Smad interactions.

Primary hepatocytes from Smad binding deficient mSnoN knock-in animals or wild-type animals were stimulated with 1 ng/ml TGFβ for the indicated times. Relative concentrations of either SnoN or Smad4 co-immunoprecipitated with Smad2/3 (A and B) as well as phosphorylated and total Smad2 (C and D) in whole cell lysates were determined by quantitative immunoblotting (data points). Trajectories (dashed lines) represent model predictions for mSnoN (grey) or wt cells (black). Standard deviations of experimental data points were calculated from triplicates (see Section 6.5.2).

explanation, i.e., very rapid degradation of constitutively expressed SnoN mRNA and thus swift termination of constitutive SnoN protein synthesis (without stimulus-enhanced degradation), seems unlikely in light of RNA experiments done by Peter Nickel under the same conditions (1 ng/ml TGFβ + actinomycin D). He found that the SnoN mRNA half-life in primary mouse hepatocytes is ~ 100 min, which is inconsistent with a half-maximal protein decay in actinomycin-treated cells within 50 min. SnoN measurements were also performed in cells where Smad signalling was terminated 50 min after TGFβ stimulation by addition of the TGFβ receptor inhibitor SB431542 (Fig. 6.8D). A very slow SnoN protein decay was observed after inhibitor treatment, most likely because Smad2-mediated SnoN degradation was blocked in addition to Smad2-mediated SnoN transcription.

The core components of TGFβ signalling such as the TGFβ receptor, Smad2 and Smad4 are all characterised by relatively stable mRNAs and/or proteins (Tables F.1 – F.3). Thus, the alterations in Smad signalling done in the presence of actinomycin D are likely to reflect at least in part transcriptional feedback regulation, since actinomycin can be assumed to mainly inhibit Smad-induced de novo transcription. Actinomycin D had little impact on Smad2 phosphorylation, but strongly affected the amount of Smad4 co-immunoprecipitated with Smad2. The observation that SnoN, a putative feedback regulator in hepatocytes, rapidly decayed in actinomycin-treated cells suggested that SnoN might be responsible for altered Smad heterotrimerisation. In Section 6.6, the impact of SnoN on Smad signalling was therefore investigated numerically and the model predictions were then confirmed experimentally.