Re: [SOSlib-devel] variable compartment bug: example[Scanned]

[Andrew F. <af...@ph...>]

Rainer
=20
> do you know of any examples of models with a compartment with=20
> spatialDimension zero or hasOnlySubstanceUnits true?
>

Darren Wilkinson created a stoctastic test suite which consists of
models
of this form.  They should operate correctly in a deterministic
simulator as well.
=20
> can a compartment of dimension zero have a size?

No

>=20
> would we currently handle such structures correctly, at least=20
> for fixed compartment sizes?

I don't see why not

yours Andrew

> -----Original Message-----
> From: Rainer Machne [mailto:ra...@tb...]=20
> Sent: 20 June 2006 15:48
> To: Andrew Finney
> Cc: SOSlib Development
> Subject: Re: [SOSlib-devel] variable compartment bug: example[Scanned]
>=20
>=20
> >>> The units of the species are of the form substance/size units
>=20
> > "The units of the species are of the form substance/size=20
> units (i.e.,=20
> > concentration units, using a broad definition of=20
> concentration) if the=20
> > compartment's spatialDimensions is non-zero and=20
> hasOnlySubstanceUnits=20
> > has the value ``false''. The units of the species are of the form=20
> > substance if spatialDimensions is zero or hasOnlySubstanceUnits has=20
> > the value ``true''. The units of substance are those defined in the=20
> > substanceUnits, and the size units are those given in the=20
> > spatialSizeUnits field."
> > ...
>=20
> oops, i guess it was some freudian process that blinded me to=20
> the continuation of this sentence :)
>=20
> do you know of any examples of models with a compartment with=20
> spatialDimension zero or hasOnlySubstanceUnits true?
>=20
> can a compartment of dimension zero have a size?
>=20
> would we currently handle such structures correctly, at least=20
> for fixed compartment sizes?
>=20
> BTW, there are several models with variable compartments in=20
> the test suite, all in the group rulesForParametersAndCompartments.
>=20
> Could it be that SOSlib only succeeded because the program=20
> used to generate tested results also does it wrong?
>=20
> Rainer
>=20
>=20
> On Tue, 20 Jun 2006, Andrew Finney wrote:
>=20
> > Rainer
> >
> >>>> of course you have to check if the species in your kinetic
> >> law have
> >>>> dimension conc or substance!
> >>>> sbml allows the usage of both, doesn't it?
> >>>
> >>> See section 4.6.4 here
> >>>
> >>=20
> http://sbml.org/specifications/sbml-level-2/version-1/html/sbml-level
> >> -
> >>> 2.html#SECTION00046000000000000000
> >>>
> >>> "
> >>> The units of the species are of the form substance/size units
> >>>
> >>> ...
> >>>
> >>> The units of the species are used in the following ways:
> >>>
> >>> ...
> >>>    * The species identifier has these units when it appears as a=20
> >>> numerical quantity in a mathematical formula expressed in MathML=20
> >>> (discussed in Section 3.6.3).
> >>> "
> >>>
> >>> So i guess, the units of species in math expressions are=20
> always in=20
> >>> concentrations!
> >>>
> >
> > WRONG
> >
> > They are either substance or concentration units depending on the=20
> > attributes of the species element.
> >
> > "The units of the species are of the form substance/size=20
> units (i.e.,=20
> > concentration units, using a broad definition of=20
> concentration) if the=20
> > compartment's spatialDimensions is non-zero and=20
> hasOnlySubstanceUnits=20
> > has the value ``false''. The units of the species are of the form=20
> > substance if spatialDimensions is zero or hasOnlySubstanceUnits has=20
> > the value ``true''. The units of substance are those defined in the=20
> > substanceUnits, and the size units are those given in the=20
> > spatialSizeUnits field."
> > ...
> >
> > "The units of the species are used in the following ways:
> >
> >
> > The species initialConcentration field has these units (see Section=20
> > 4.6.3).
> >
> > The species identifier has these units when it appears as a=20
> numerical=20
> > quantity in a mathematical formula expressed in MathML=20
> (discussed in=20
> > Section 3.6.3).
> >
> > The math field of an AssignmentRule structure determining=20
> the species'
> > quantity (see Section 4.8.2) has these units.
> >
> > In RateRule structures that set the rate of change of the species'
> > quantity (Section 4.8.3), the units on the rule's math=20
> field are the=20
> > units of the species divided by the built-in time units. "
> >
> > yours Andrew
> >
> >
> >
> >
> >> -----Original Message-----
> >> From: sbm...@li...
> >> [mailto:sbm...@li...] On=20
> Behalf Of=20
> >> Rainer Machne
> >> Sent: 20 June 2006 15:26
> >> To: Stefan Mueller
> >> Cc: sbm...@li...
> >> Subject: Re: [SOSlib-devel] variable compartment bug:=20
> >> example[Scanned]
> >>
> >>
> >>
> >> This
> >>
> >>> current ODE: d[S1]/dt =3D f(S1, S2, ...) / V planned ODE: d=20
> S1 /dt =3D=20
> >>> f(S1, S2, ...)
> >>
> >> should of course be
> >>
> >> current ODE: d[S1]/dt =3D f([S1], [S2], ...) / V planned ODE: d
> >> S1 /dt =3D f([S1], [S2], ...)
> >>
> >>
> >> On Tue, 20 Jun 2006, Rainer Machne wrote:
> >>
> >>> stefan
> >>>
> >>>
> >>>>> divide k by V(0) and the whole kinetic law again by V^2 for 2nd=20
> >>>>> order reactions
> >>>>
> >>>> yes.
> >>>> and that's exactly the same what you are planning, if you
> >> want to use
> >>>> conc ans subst in parallel.
> >>>> only the point of conversion differs.
> >>>
> >>> not really. i wasn't planning any conversion of the math,
> >> except for:
> >>>
> >>> current ODE: d[S1]/dt =3D f(S1, S2, ...) / V planned ODE: d=20
> S1 /dt =3D=20
> >>> f(S1, S2, ...)
> >>>
> >>> then we have two data->value arrays, let's say
> >>>
> >>> data->odeValues =3D S1, S2, S3, etc.
> >>> and
> >>> data->concValues =3D [S1], [S2], [S3] etc.
> >>>
> >>>
> >>> The first array `odeValues' corresponds to the SUNDIALS
> >> values in NVector y.
> >>> The second array `concValues' will be used in evaluateAST, for=20
> >>> evaluation the math.
> >>>
> >>>
> >>> This would be equivalent of course, with just replacing=20
> all [S1] in=20
> >>> all math expressions with S1/volume, and only using one array
> >>> data->odeValues
> >>>
> >>> While the latter approach might be easier to implement it has a=20
> >>> probably the not so small drawback, that for each species
> >> in each math
> >>> expression the evaluation needs one more operation, i.e.
> >> the division,
> >>> while this division could be done only once for each
> >> species whenever
> >>> data->odeValues is updated (i.e. e.g. in the
> >> SOSlib/SUNDIALS funciton f).
> >>>
> >>> However, this might still only be correct, if compartments are=20
> >>> determined by assignment rules.
> >>> If they are determined by rate rules, then the replacement
> >> S1 -> S1/V
> >>> would also be required for constructing the jacobi matrix.
> >>>
> >>>> of course you have to check if the species in your kinetic
> >> law have
> >>>> dimension conc or substance!
> >>>> sbml allows the usage of both, doesn't it?
> >>>
> >>> See section 4.6.4 here
> >>>
> >>=20
> http://sbml.org/specifications/sbml-level-2/version-1/html/sbml-level
> >> -
> >>> 2.html#SECTION00046000000000000000
> >>>
> >>> "
> >>> The units of the species are of the form substance/size units
> >>>
> >>> ...
> >>>
> >>> The units of the species are used in the following ways:
> >>>
> >>> ...
> >>>    * The species identifier has these units when it appears as a=20
> >>> numerical quantity in a mathematical formula expressed in MathML=20
> >>> (discussed in Section 3.6.3).
> >>> "
> >>>
> >>> So i guess, the units of species in math expressions are=20
> always in=20
> >>> concentrations!
> >>>
> >>> Rainer
> >>>
> >>>
> >>>
> >>>
> >>>
> >>>>
> >>>> i mean, every species is in a well-defined compartment, isn't it?
> >>>> so you can use the formula [S]=3DS/V whenever you want...
> >>>>
> >>>> i would eliminate conc as early as possile, but that's a=20
> matter of
> >>> taste.
> >>>
> >>>
> >>>
> >>> On Tue, 20 Jun 2006, Stefan Mueller wrote:
> >>>
> >>>> hi rainer!
> >>>>
> >>>>> stefan
> >>>>>
> >>>>>>> Currently the ODE would be:
> >>>>>>>
> >>>>>>> d[S1]/dt =3D k1 * [S1] / V1
> >>>>>>
> >>>>>> this is not valid for variable compartments.
> >>>>>
> >>>>> i said that it's invalid for variable compartments some
> >> lines below:
> >>>>>>> which is only valid for constant compartments V1 and V2.
> >>>>>>
> >>>>>> please don't use it for comparison!
> >>>>>
> >>>>> it's not used for comparison, but it is the current
> >> implementation
> >>>>> in SOSlib!
> >>>>
> >>>> i said this just to avoid any additional misunderstanding.
> >>>> the whole topic seems to be quite confusing...
> >>>>
> >>>>>>> At least for the second ODE, depending in two different
> >> volumes,
> >>>>>>> the conversion to ODEs depending on substance instead of=20
> >>>>>>> concentration is not as simple as
> >>>>>>>
> >>>>>>>> dS1/dt =3D k * [S1] * [S2] / V(0) * V =3D k/V(0) * S1 * S2 / =
V
> >>>>>>
> >>>>>> why dou you refer to this law for a reaction of order 2?
> >>>>>
> >>>>> sorry, at first i wanted to make the second reaction in=20
> V2 of 2nd=20
> >>>>> order, e.g.,
> >>>>>
> >>>>> in V2:
> >>>>> S2 + S3 -> P: k2*[S2]*[S3]
> >>>>>
> >>>>>>> as your example before,
> >>>>>>> but for S2 I get:
> >>>>>>>
> >>>>>>> dS2/dt =3D k1/V1(0) * S1/V1 * V2  -  k2/V2(0) * S2
> >>>>>>>
> >>>>>>> Right?
> >>>>>>
> >>>>>> no!
> >>>>>> dS2/dt =3D k1/V1(0) * S1 - k2/V2(0) * S2 it's as simple as =
that!
> >>>>>> (since dS2/dt =3D -dS1/dt as far as reaction S1->S2 is=20
> concerned.)=20
> >>>>>> that's the advantage of kinetic laws (and consequently=20
> odes) for=20
> >>>>>> substances!
> >>>>>
> >>>>> of course! sorry, i got quite confused obviously.
> >>>>>
> >>>>> and for the above second order reaction in V2 it would be:
> >>>>>
> >>>>>   dS2/dt =3D k1/V1(0) * S1 - k2/V2(0) * S1 * S3 / V2
> >>>>>
> >>>>> Right this time?
> >>>>
> >>>> yes :)
> >>>>
> >>>>>>> Thus I would really have to check the type of reaction
> >> and edit a
> >>>>>>> lot of kinetic laws before constructing the ODEs.
> >>>>>>
> >>>>>> the above example shows that this is not necessary.
> >>>>>> maybe there are other examples...
> >>>>>
> >>>>> so you say, we could convert kinetic laws that don't explicitly=20
> >>>>> contain the volume by simply differentiating these=20
> three cases and
> >>>>>
> >>>>>>>> order 1:
> >>>>>>>> S1->junk: k*[S1]
> >>>>>>>> dS1/dt =3D k * [S1] / V(0) * V =3D k/V(0) * S1
> >>>>>
> >>>>> divide k by V(0) for 1st order reactions
> >>>>>
> >>>>>>>> order 2:
> >>>>>>>> S1+S2->bla: k*[S1]*[S2]
> >>>>>>>> dS1/dt =3D k * [S1] * [S2] / V(0) * V =3D k/V(0) * S1 * S2 / =
V
> >>>>>
> >>>>> divide k by V(0) and the whole kinetic law again by V^2 for 2nd=20
> >>>>> order reactions
> >>>>>
> >>>>>>>> order 3:
> >>>>>>>> S1+S2+S3->wow: k*[S1]*[S2]*[S3]
> >>>>>>>> dS1/dt =3D k * [S1] * [S2] * [S3] / V(0) * V =3D k/V(0) *
> >> S1 * S2 *
> >>>>>>>> S3 /
> >>>>>
> >>>>> V^2
> >>>>>
> >>>>> divide k by V(0) and the whole kinetic law again by V^2 for 2nd=20
> >>>>> order reactions
> >>>>
> >>>> yes.
> >>>> and that's exactly the same what you are planning, if you
> >> want to use
> >>>> conc ans subst in parallel.
> >>>> only the point of conversion differs.
> >>>>
> >>>> i mean, every species is in a well-defined compartment, isn't it?
> >>>> so you can use the formula [S]=3DS/V whenever you want...
> >>>>
> >>>> i would eliminate conc as early as possile, but that's a
> >> matter of taste.
> >>>>
> >>>>> but there are also enzymatic reactions, with in often
> >> quite complex
> >>>>> kinetic laws. we would really have to do unit checking,
> >> recognizing
> >>>>> the implemented reaction mechanisms etc.
> >>>>
> >>>> of course you have to check if the species in your kinetic
> >> law have
> >>>> dimension conc or substance!
> >>>> sbml allows the usage of both, doesn't it?
> >>>>
> >>>>> also, we don't really know in which compartment=20
> reactions happen,=20
> >>>>> but could just try to get this information from the
> >> reactants' compartments.
> >>>>
> >>>> this does not matter, for 2 reasons:
> >>>>
> >>>> 1. if the kinetic law is specified correctly (including the=20
> >>>> compartment size), e.g. A+B->C: k*V*[A]*[B], the the
> >> compartment size
> >>>> is considered automatically.
> >>>>
> >>>> 2. the conversion [S]=3DS/V does not depend on reactions,
> >> but only on
> >>>> species (and their compartments).
> >>>>
> >>>>> so i guess, it would be much easier to
> >>>>>
> >>>>> a) use ODEs for substances but dependent on
> >> concentrations, here the
> >>>>> only thing we need to do is to leave out the compartment
> >> division at
> >>>>> the end of Species_odeFromReactions and have a double
> >> bookkeeping of
> >>>>> amounts and concentrations
> >>>>
> >>>> to my taste, double book keeping is dangerous and produces
> >> overhead.
> >>>>
> >>>>> and
> >>>>>
> >>>>> b) just ignore WRONG models that have a variable=20
> compartment but=20
> >>>>> don't have the compartment explicitly in their kinetic laws (or=20
> >>>>> check for such cases and issue an error);
> >>>>>
> >>>>> isn't it an issue of correct model writing rather then
> >> correct model
> >>>>> solving?
> >>>>
> >>>> i agree.
> >>>>
> >>>>> However, another problem might be that to calculate in amounts=20
> >>>>> rather then in concentrations also changes the use of absolute=20
> >>>>> errors and maybe also the general the performance of the
> >> integrator, i guess. ??
> >>>>
> >>>> abs errors might change, but that's not an issue.
> >>>> rel errors should stay the same.
> >>>> (in one-compartment models everything is just multiplied by the
> >>>> volume.)
> >>>>
> >>>>> maybe we will have to really implement both in parallel,
> >>>>
> >>>> i hope not :)
> >>>>
> >>>>> * ODEs for concentrations/time if no variable=20
> compartments are in=20
> >>>>> the model
> >>>>> * ODEs for ammounts/time if any compartment is variable
> >>>>>
> >>>>> ?
> >>>>>
> >>>>> Rainer
> >>>>
> >>>> maybe we discuss the next round via telephone :)
> >>>>
> >>>> cheers, s
> >>>>
> >>>>
> >>>>>>> Thus I think it just doesn't work with kinetic laws=20
> that don't=20
> >>>>>>> include the compartment explicitly!
> >>>>>>
> >>>>>> i don't know exactly what you mean.
> >>>>>> of course kinetic laws with explicit compartment are
> >> "nicer", but
> >>>>>> you can always transform any kinetic law to a "nice" one...
> >>>>>>
> >>>>>>> But if they include the compartment we can write ODEs for=20
> >>>>>>> substances but have concentrations in the right hand
> >> side of the
> >>>>>>> ODEs. We then just need to a second array for the
> >> concentrations, as discussed before.
> >>>>>>
> >>>>>> if we still use concentrations is a separate problem,
> >> independent
> >>>>>> of which we discussed until now.
> >>>>>> it's a matter of taste, if we eliminate the=20
> concentrations  from=20
> >>>>>> the kinetic laws and get a (possibly nonlinear)
> >> dependence on the
> >>>>>> compartments OR if we keep concentrations and always=20
> have linear=20
> >>>>>> dependence on the comp.
> >>>>>>
> >>>>>>> Rainer
> >>>>>>
> >>>>>> cheers,
> >>>>>> stefan.
> >>>>>>
> >>>>>>> On Mon, 19 Jun 2006, Stefan Mueller wrote:
> >>>>>>>> let me summarize:
> >>>>>>>> kinetic laws have units substance/time, and they depend on=20
> >>>>>>>> concentrations or substances.
> >>>>>>>>
> >>>>>>>> we SHOULD write down odes for substances, and we COULD
> >> eliminate
> >>>>>>>> concs completely.
> >>>>>>>> (of course, we would have to convert initial concs to
> >>>>>>>> substances.)
> >>>>>>>>
> >>>>>>>> in any case, odes for substances depend on the
> >> compartment size
> >>>>>>>> (linearly or according to the order of the reaction):
> >>>>>>>>
> >>>>>>>> order 1:
> >>>>>>>> S1->junk: k*[S1]
> >>>>>>>> dS1/dt =3D k * [S1] / V(0) * V =3D k/V(0) * S1
> >>>>>>>>
> >>>>>>>> order 2:
> >>>>>>>> S1+S2->bla: k*[S1]*[S2]
> >>>>>>>> dS1/dt =3D k * [S1] * [S2] / V(0) * V =3D k/V(0) * S1 * S2 / =
V
> >>>>>>>>
> >>>>>>>> order 3:
> >>>>>>>> S1+S2+S3->wow: k*[S1]*[S2]*[S3]
> >>>>>>>> dS1/dt =3D k * [S1] * [S2] * [S3] / V(0) * V =3D k/V(0) *
> >> S1 * S2 *
> >>>>>>>> S3 /
> >>>>>>>> V^2
> >>>>>>>>
> >>>>>>>> ...
> >>>>>>>>
> >>>>>>>> if the kinetic laws included the compartment
> >> explicitly (which is
> >>>>>>>> desirable), we would have the following situation:
> >>>>>>>>
> >>>>>>>> order 1:
> >>>>>>>> S1->junk: r*[S1]*V
> >>>>>>>> dS1/dt =3D  r * [S1] * V=3D r * S1
> >>>>>>>>
> >>>>>>>> order 2:
> >>>>>>>> S1+S2->bla: r*[S1]*[S2]*V
> >>>>>>>> dS1/dt =3D r * [S1] * [S2] * V=3D r * S1 * S2 / V
> >>>>>>>>
> >>>>>>>> order 3:
> >>>>>>>> S1+S2+S3->wow: r*[S1]*[S2]*[S3]*V
> >>>>>>>> dS1/dt =3D r * [S1] * [S2] * [S3] * V =3D r * S1 * S2 * S3 / =
V^2
> >>>>>>>>
> >>>>>>>> as stated earlier, odes for substances have got another=20
> >>>>>>>> advantage. for the compartment size we could use=20
> either a rate=20
> >>>>>>>> rule or an assignment
> >>>>>>>> rule: dV/dt =3D f(t, V, ...)
> >>>>>>>> or
> >>>>>>>> V=3D f(t, V, ...)
> >>>>>>>>
> >>>>>>>> cheers,
> >>>>>>>> stefan.
> >>>>>>>>
> >>>>>>>> Am Montag 19 Juni 2006 11:17 schrieb Stefan Mueller:
> >>>>>>>>> andrew:
> >>>>>>>>>
> >>>>>>>>> for the SBML system
> >>>>>>>>> 	S1 in C
> >>>>>>>>> 	S1->junk : k * [S1]
> >>>>>>>>> 	dC/dt =3D 1
> >>>>>>>>>
> >>>>>>>>> the right math is:
> >>>>>>>>> 	dS1/dt =3D k * [S1] / C(0) * C
> >>>>>>>>> 	dC/dt =3D 1
> >>>>>>>>> with initial cond:
> >>>>>>>>> 	S1(0) =3D [S1](0) * C(0)
> >>>>>>>>> 	C(0)
> >>>>>>>>>
> >>>>>>>>> notation:
> >>>>>>>>> 	S1 ... substance
> >>>>>>>>> 	[S1] ... conc
> >>>>>>>>>
> >>>>>>>>> the "proposed solution" (who proposed this?) has the right=20
> >>>>>>>>> units, but is wrong otherwise... :)
> >>>>>>>>>
> >>>>>>>>> as rainer and me already discussed, everything=20
> would be much=20
> >>>>>>>>> clearer, if the kinetic law included the compartment
> >> explicitly:
> >>>>>>>>> 	S1 in C
> >>>>>>>>> 	S1->junk : C * r * [S1]
> >>>>>>>>> 	dC/dt =3D 1
> >>>>>>>>>
> >>>>>>>>> with:
> >>>>>>>>> 	r ... the "real" chemical rate constant.
> >>>>>>>>>
> >>>>>>>>> rainer:
> >>>>>>>>>
> >>>>>>>>> chemists always knew that rate constants are not
> >> constant with
> >>>>>>>>> respect to pH, T, ... as a physicist i have to be fair! :)
> >>>>>>>>>
> >>>>>>>>> to conclude: the math is not really tricky.
> >>>>>>>>> odes for substances are better suited to handle=20
> both multiple=20
> >>>>>>>>> and variable compartments. of course, we have to=20
> take care of=20
> >>>>>>>>> many implementation details...
> >>>>>>>>>
> >>>>>>>>> cheers,
> >>>>>>>>> stefan.
> >>>>>>>>
> >>>>>>>> _______________________________________________
> >>>>>>>> sbmlsolver-devel mailing list
> >>>>>>>> sbm...@li...
> >>>>>>>> https://lists.sourceforge.net/lists/listinfo/sbmlsolver-devel
> >>>>>>>
> >>>>>>> _______________________________________________
> >>>>>>> sbmlsolver-devel mailing list
> >>>>>>> sbm...@li...
> >>>>>>> https://lists.sourceforge.net/lists/listinfo/sbmlsolver-devel
> >>>>>>
> >>>>>> _______________________________________________
> >>>>>> sbmlsolver-devel mailing list
> >>>>>> sbm...@li...
> >>>>>> https://lists.sourceforge.net/lists/listinfo/sbmlsolver-devel
> >>>>>
> >>>>> _______________________________________________
> >>>>> sbmlsolver-devel mailing list
> >>>>> sbm...@li...
> >>>>> https://lists.sourceforge.net/lists/listinfo/sbmlsolver-devel
> >>>>
> >>>>
> >>>> _______________________________________________
> >>>> sbmlsolver-devel mailing list
> >>>> sbm...@li...
> >>>> https://lists.sourceforge.net/lists/listinfo/sbmlsolver-devel
> >>>>
> >>>
> >>
> >>
> >> _______________________________________________
> >> sbmlsolver-devel mailing list
> >> sbm...@li...
> >> https://lists.sourceforge.net/lists/listinfo/sbmlsolver-devel
> >>
> >
> >
> > _______________________________________________
> > sbmlsolver-devel mailing list
> > sbm...@li...
> > https://lists.sourceforge.net/lists/listinfo/sbmlsolver-devel
> >
>=20