m = pwAddY(m, *ID, *rhs, errorModelRhs, noiseType, unit, name,
description, alternativeIDs, designerProps, classname)
m = pwAddY(m, interactive)
Adds an observable to a model. An observable connects a model to an
NEW SYNTAX SINCE PW 3.0:
pwAddY(m, rhs, ID, scalingParameter, errorModelRhs, noiseType, ...)
pwAddY(m, ID, rhs, errorModelRhs, noiseType, ...)
Arguments for pwAddY
ID ID of the observable. Should be different
than any other ID, e.g. of a dynamic variable X.
rhs Right-hand side of the observation
May include x, k, a, and s IDs.
errorModelRhs Default: 'y * 0.1 + max(y) * 0.05'
noiseType Currently only 'Gaussian' supported.
name descriptive name
description Description of the observation
alternativeIDs Cell array of alternative IDs which will be used for
automatic mapping when adding external data.
designerProps for internal use
Alternative function call for interactive specification:
interactive If true, the user is prompted to enter the properties
into a user interface.
Observables Y are functions of dynamical variables X, dynamical parameters K,
algebraic assignments A, and scaling factors S. They correspond to measurable
quantities in an experiment. The right-hand side of the observation function
is specified as a mathematical expression of the IDs belonging X, K, and S.
In PW 2.* an implicit scaling parameter was used for each observable.
This is not the case any more since PW 3.0.
Now, a scaling between rhs and measurements has to be compensated by
an S parameter introduced into the rhs directly, e.g.
pwAddY(m, 'pS_obs', 'scale_pS_obs * (pS + 2 * pStat_pS)');
pwAddS(m, 'scale_pS_obs', 1);
This notation is more transparent and explicit.
Each data point will be weighted during model fitting depending on
its error bar, i.e. standard deviation. The standard deviation is
either calculated using an error model, estimated on a spline-based
approximation or given explicitly by the user in the data file (see pwAddData).
A default error model is always required, e.g. in order to simulate
synthetic data. To specify Gaussian noise with a standard deviation
of 10% relative to y plus 5% absolute (relative to max(y) over all y):
m = pwAddY(m, 'pS_obs', 'scale_pS_obs * (pS + 2 * pS_pS)', '0.1*y + 0.05*max(y)');
When using max(yAllStimuli)instead of max(y), the maximum of all
y values of all stimuli will be calculated.
Fitting with logarithmic data
In Kreutz et al. 2007 it has been suggested to use log-normally distributed errors
e.g. in Western blotting. Not taking into account systematic effects e.g.
of time or gel lane, the measurement model reads:
yMeas = beta0 + beta1 * yTrue^beta2 * eta with eta ~ e^N(0,sigma)
We assume that the measured data have no offset, i.e. beta0 = 0.
Taking the logarithm yields:
log(yMeas) = log(beta1) + beta2 * log(yTrue) + log(eta)
beta2 is usually 1 (linear dependency between true and measured y) leading to:
log(yMeas) = log(beta1) + log(yTrue) + log(eta).
Note that the scaling factor beta1 is now an additive offset parameter.
Since log(eta) ~ N(0, sigma), for log(yMeas) a normal error model can be
assumed, distributed as N(0, sigma).
Therefore, if you take the logarithm of an experimentally measured
observable A, the following observable has to be added to the model
(assuming log-normally distributed errors with mu = 0 and sigma = 0.2):
m = pwAddY(m, 'logStat_obs', 'log(A) + log(beta1)', '0.2');
Simultaneously, we want to fit locally the scaling parameter which
enters the right-hand side as the additive constant beta1. With an initial
guess of 1 we use:
pwAddS(m, 'beta1', 1, 'local');
If you prefer log10, just use 1.0 / log(10.0) * log(A)
or 0.434294481903252 * log(A).
This corresponds to a changed basis (=log10(A)).
Note that this also affects sigma.
Please compare with
C. Kreutz, M.M. Bartholome Rodriguez, T. Maiwald, M. Seidl, H.E. Blum, L. Mohr, J. Timmer.
An error model for protein quantification. Bioinformatics 23, 2007, 2747-2753