sbmlsolver-devel

Re: [SOSlib-devel] comments on SBML L2 V2

[Rainer M. <ra...@tb...>]

Hi Michael,

This email also goes to the SOSlib developers' email list.

On Mon, 4 Sep 2006, Michael Hucka wrote:

> Andrew recently wrote the following to me in a separate
> message about the L2v2 draft specification:
>
>  afinney> the Example of the translation from
>  afinney> multicompartment models to odes:
>  afinney>
>  afinney> This model is correct as it stands but this type
>  afinney> of translation doesn't work when the compartment
>  afinney> volume varies i.e. is defined by an ODE.  In that
>  afinney> type of system the species odes need to be in
>  afinney> substance/time units.  The species concentrations
>  afinney> can then be computed by simple 'rules' applied to
>  afinney> intermediate species variables in substance
>  afinney> units.  The SOSlib team ran into this issue a
>  afinney> month or 2 ago.  I can't explain the math but we
>  afinney> could probably get them to explain this if you
>  afinney> want real proof.
>  afinney>
>  afinney> Perhaps we can have a variable compartment model
>  afinney> translated into Odes in addition to the existing
>  afinney> model
>  afinney>
>  afinney> either way there needs to some text describing
>  afinney> the limitations of the approach that's used in
>  afinney> the current text.
>
> I think the model in question in the examples section was
> contributed by Nicolas.  I'm afraid I'll run out of time
> to make all the remaining changes to the spec, so could I
> ask for people's help in doing as Andrew suggests above?

The example in section 5.6 has a variable compartment V1:
http://sbml.org/specifications/sbml-level-2/version-1/html/sbml-level-2.html#SECTION00056000000000000000

So I think the given translation to ODEs is wrong in the current state!! 
Or am I wrong? This topic confused me back then. Maybe Stefan or Luki can 
check section 5.6 and the following short explanation:

1) First the constant K3 (and all other Ks) is already corrected for some 
initial volume to convert the kinetic law to substance/time.

K3 = K3original * V1(0)
where K3original was in units [1/time] and K3 is in [volume/time], so 
K3*S3, with S3 as a concentration [substance/volume], yields the required 
[substance/time].

So K3 is actually a volume-dependent variable!!
The equation should rather have an explicit volume term, i.e.

K3original * V1 * S1

The volume corrections in kinetic laws should in general be done 
explicitly and it might be a good thing to note this in the specs 
compartment section.

2) Second the last ODE dS3/dt is wrong. A correct math would require a 
d[S3]/dV term to account for the influence of volume change on the 
concentration.

Imagine an initial concentration of 1 for S3(0) and a Volume V1 that has 
doubled at time t, V1(t) = 2 * V1(0). The ODE given would just add the 
changes by reactions but not changes to this initial condition by the 
volume change which at time t would only be 0.5. A term dS3/dV would 
account for this.

Fortunately there was an easier solution which Andrew mentioned above and 
which is outlined shortly at 
http://www.tbi.univie.ac.at/wiki/index.php/SOSlibDevel#Variable_Compartment_Sizes

1) write ODEs in substance/time, i.e. don't divide by the compartment, the 
left hand side of the ODEs are dXamount/dt rather then dX/dt; the ODE 
solver will the calculate Xamount.
2) write assignment rules for species concentrations which are required at 
the right hand side of ODEs:
X = Xamount / Volume.

Let me know if we can contribute something specific to the new draft.

Rainer



>
> The possible outcomes are as follows, in decreasing value of
> goodness:
>
> 1. Someone writes a 2nd model (and corresponding
>   explanation) that handles the variable volume case.  We
>   add this model as an additional example in the examples
>   section, and also modifies the description of the
>   existing constant-volume model to explain it's a common
>   but special case, it's limited, doesn't handle varying
>   volumes, etc.
>
> 2. Someone modifies the existing model to make it handle
>   varying volumes, and we leave just the one model, because
>   it's the most general case.
>
> 3. We decide we don't have time to create a new general
>   example, and instead only edit the existing
>   constant-volume model to explain its limitations.
>
> 4. I add a sentence to the existing constant-volume model
>   example description mentioning that it's only for
>   constant volumes "and this does not properly handle the
>   more general case".
>
> Let me know if anyone can do something.
>
> MH
>

[SOSlib-devel] Request for help with a multicompartment problem

[Michael H. <mh...@ca...>]

Rainer and the rest,

Thanks for finding those errors in the examples in the SBML
spec; we'll try to get them out.

But since you have been thinking about this topic, I'd like
to ask you how you think the following class of problems can
be solved.  It has been confounding Nicolas and myself.

Suppose we have a reaction system consisting of 3
compartments, with species A in Va, B in Vb, and C in Vc
(I'm using the volumes as the names of the compartments),
and the reaction is

   A + B -> C

The compartments have constant volume; no varying volumes
here.  Suppose further that the reaction is said to obey the
simplest possible kind of kinetics, a kind of mass action
(but a weird one because multiple compartments are
involved).  So let's say the rate is

  k * [A] * [B]

Here, k must have units of a 2nd order rate constant, which
is 1/(Molar * time) = (volume)/(substance*time).  As initial
values, we might be given the values of concentrations [A],
[B], [C], or we might be given the substance amounts; it
won't matter as we'll see below.  The key point is that the
reactants are in different compartments.

Now we want to construct the ODEs.

First let's start naively just to make sure all the grounds
are covered.  The first attempt we might make is to say that
the rates of change are like in traditional kinetics:

  d[A]/dt = -k * [A] * [B]
  d[B]/dt = -k * [A] * [B]
  d[C]/dt =  k * [A] * [B]

But this can easily be shown to be wrong.  Borrowing an
example that Nicolas has written for the l2v2 spec, let's
say that volume Vc = 2*Va (i.e., Vc is twice the volume of
Va).  Then if we calculate the change of substance amounts
(moles or items) in compartments Va and Vc (and here I'm
using "A" to stand for the amount of substance of species A,
"C" to stand for the amount of substance of species C):

  dA/dt = Va*d[A]/dt = Va * (-k * [A] * [B]^2)

  dC/dt = Vc*d[C]/dt = Vc * (k * [A] * [B]^2)
                     = 2 * Va * (k * [A] * [B]^2)

So we just created matter out of nothing.  This simply
proves what we already knew, which is that we have to
express the rate in terms of substance changes, not
concentration changes.

OK, fine, so let's try to do that.  Since we know the
compartment volumes and we can know the concentrations of
each species, we have the relationships

  [A] = A/Va
  [B] = B/Vb
  [C] = C/Vc

We then substitute this into the ODEs:

  d[A]/dt = -k * (A/Va) * (B/Vb)
  d[B]/dt = -k * (A/Va) * (B/Vb)
  d[C]/dt =  k * (A/Va) * (B/Vb)

So far this is perfectly valid, but it didn't help anything;
the overall expressions are still in terms of
concentrations.  We need to make them to be about species
amount, not concentrations.

But how?

Attempt #1: Multiply each ODE by the volume in which that
species is located.  Let's take for example d[C]/dt, and
multiply it by the volume in which [C] is located:

  dC/dt = Vc*d[C]/dt = k * (A/Va) * (B/Vb) * Vc

This seems reasonable; after all, C is located in Vc, and
Vc*[C] must give me number of molecules of C.  But it can't
be right.  Suppose I change the model by doubling the size
of Vc, without changing anything else in the system.
According to the above, the rate of appearance of molecules
of C will double if I do that, which is impossible because
actually the number of molecules must be solely determined
by the reactants.

Attempt #2: The conclusion from the above must be that
multiplying by Vc is the wrong way around.  The core rate
expression must be reformulated in terms of substance per
time, not concentration per time, and *then* concentrations
can be made out of it by dividing by the volumes as needed.
So let's try to do that.

  new rate expr = k * ([A]*Va) * ([B]*Vb)

This gets rid of the concentration terms, but now the units
are wrong, because they are (volume * substance)/time rather
than substance/time.  Hmmm.  Could this be a clue?  Maybe we
just need to divide by a volume?  That seems promising.  So:

  dA/dt = ( -k * ([A]*Va) * ([B]*Vb) )/Va
  dB/dt = ( -k * ([A]*Va) * ([B]*Vb) )/Vb
  dC/dt = (  k * ([A]*Va) * ([B]*Vb) )/Vc

The units look right.  But wait, this has the same problem
as before: if Vc changes and nothing else in the system
changes, the rate of appearance of C will change.

Attempt #3: Maybe we simply divided by the wrong volume in
the equations above; maybe it should be always Va or Vb.
Oh, but, there are two choices: Va and Vb.  If we can only
use one, how do we pick which one we use?

Attempt #4: Ok, so the volume division idea isn't working.
What if the issue is not about dividing by volumes, but
instead, the problem is that the units of k need to be
changed?  Maybe what we really need is

  dA/dt = -k' * ([A]*Va) * ([B]*Vb)
  dB/dt = -k' * ([A]*Va) * ([B]*Vb)
  dC/dt =  k' * ([A]*Va) * ([B]*Vb)

Then k' can be used to sort out the units.  Except.... what
*should* be the correction applied to k?  How can we create
k' systematically?  What if the rate law is something
entirely different -- without a k', or with other constants?
What will the correction be then?


So at this point, we're stuck.  The conclusion must be one
of the following:

  Option (a): Conclude that attempt #4 is the right way, and
    we need to work out how to construct a proper correction
    term for the general case.

  Option (b): Conclude that there is a fundamental error in
    one or more of the attempts #1-4, or there is another
    way out.

  Option (c): Conclude that we must declare the situation
    impossible and require that all reactants must be in the
    same compartment.

Unfortunately, Nicolas points out that (c) is violated by
some existing models -- they already use multiple
compartments for reactants -- so people must believe it's
doable.  Option (c) is not a way out.

What is your opinion?  What is the solution?  How should
this multicompartment model be solved?  

MH

Re: [SOSlib-devel] Request for help with a multicompartment problem

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

Folks

Assuming the issue is just how to create an encoding in SBML for
multicomparment models
and not how to map from SBML to ODEs...

I'd just like to point out one obvious point...

an equation like 

d[C]/dt = k * [A] * [B]

describes a phenomena that only really works in a single compartment.
i.e. when two species are in the same compartment and therefore can freely
react with each
other the rate of reaction is proportional to the concentrations of the
species.

if A B and C are not in the same compartment the equation doesn't make any
physical sense to describe the interactions

in that case any chemical equation A + B -> C
has be decomposed into its constituent transport and conversation components
before you can move forward,

As far as SBML is concerned there 2 ways to go either you formulate a
complex rate equation which enables you to fold in the transport and
conversion components together or
you decompose the SBML reaction into simpler ones

Must go to work now

Bye

Andrew

Re: [SOSlib-devel] Request for help with a multicompartment problem

[Nicolas Le N. <le...@eb...>]

On Mon, 11 Sep 2006, Andrew Finney wrote:

> Assuming the issue is just how to create an encoding in SBML for
> multicomparment models
> and not how to map from SBML to ODEs...
>
> I'd just like to point out one obvious point...
>
> an equation like
>
> d[C]/dt = k * [A] * [B]
>
> describes a phenomena that only really works in a single compartment.
> i.e. when two species are in the same compartment and therefore can freely
> react with each
> other the rate of reaction is proportional to the concentrations of the
> species.
>
> if A B and C are not in the same compartment the equation doesn't make any
> physical sense to describe the interactions
>
> in that case any chemical equation A + B -> C
> has be decomposed into its constituent transport and conversation components
> before you can move forward,

But that is the point. Developing a model implies a balance of
abstraction Vs. biological exactitude. I recently submitted a paper
describing a model with 159 reactions. Two-third of them involve
various states of just one molecule. On the contrary, the remaining
reactions are mainly phenomenological descriptions of all the others
participants. And that's fine. Because we don't care much about
them. We just want them to behave as we think they are in the cell.

If that suits me, I can perfectly model degradation of mRNA using a
RNAse in the cytosol and the mRNA in the nucleus, without modeling the
transport of RNA. No referee will reject my model, except if it deals
specifically with RNA degration (funny enough regarding the model
cited above, we braced ourselves against the referee's criticisms on
the way we modeled second-knives and prepared alternative version of
the model. But none of the three referees - who did a very thorough
job and required a huge number of modifications - even mentioned the
issus).

IMHO, it is not our job, as SBML developers, to decide how a modeller
should formulate the model.


-- 
Nicolas LE NOVERE,  Computational Neurobiology,
EMBL-EBI, Wellcome-Trust Genome Campus, Hinxton, Cambridge, CB10 1SD, UK
Tel: +44(0)1223 494 521,  Fax: +44(0)1223 494 468,  Mob: +33(0)689218676
http://www.ebi.ac.uk/~lenov, AIM: nlenovere,  MSN: nle...@ho...

Re: [SOSlib-devel] Request for help with a multicompartment problem

[Michael H. <mh...@ca...>]

Hi Andrew,

  afinney> describes a phenomena that only really works in a
  afinney> single compartment.  i.e. when two species are in
  afinney> the same compartment and therefore can freely
  afinney> react with each other the rate of reaction is
  afinney> proportional to the concentrations of the
  afinney> species.

I suspect Nicolas will say something about it too, but I
want to mention I got an email from Bob Phair on the same
question, and he entirely contradicted the claim that this
is an unusual situation.  In his words, "The situation you
describe actually occurs frequently. I don't see anything
weird here", and goes on to give an example of extracellular
signal molecule binding to membrane receptors where you have
to deal with separate volumes.

(His solution to the problem was change the rate constant.
I'll ask if I can forward his reply around.)

What may be a bit confusing is the rate law I used,
mass-action, has assumptions about molecular interactions.
If you remove that, perhaps it won't require thinking that
the reactants have to be in the same space.

More generally, I think you would agree that modelers should
be able to express a compound rate expression that involves
whatever compartmental organization they feel is appropriate
for their problem.  (And indeed, Nicolas will attest to the
fact that models in the wild do this.)

MH

Re: [SOSlib-devel] Request for help with a multicompartment problem

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

Mike

you wrote:

> More generally, I think you would agree that modelers should
> be able to express a compound rate expression that involves
> whatever compartmental organization they feel is appropriate
> for their problem.  (And indeed, Nicolas will attest to the
> fact that models in the wild do this.)
> 

I agree.  But we don't need to get stressed if they screw up the units!

yours
Andrew Finney  

> -----Original Message-----
> From: Michael Hucka [mailto:mh...@ca...] 
> Sent: 11 September 2006 08:26
> To: Andrew Finney
> Cc: 'Rainer Machne'; 'Nicolas Le Novere'; 'SOSlib Development'
> Subject: RE: Request for help with a multicompartment problem
> 
> Hi Andrew,
> 
>   afinney> describes a phenomena that only really works in a
>   afinney> single compartment.  i.e. when two species are in
>   afinney> the same compartment and therefore can freely
>   afinney> react with each other the rate of reaction is
>   afinney> proportional to the concentrations of the
>   afinney> species.
> 
> I suspect Nicolas will say something about it too, but I
> want to mention I got an email from Bob Phair on the same
> question, and he entirely contradicted the claim that this
> is an unusual situation.  In his words, "The situation you
> describe actually occurs frequently. I don't see anything
> weird here", and goes on to give an example of extracellular
> signal molecule binding to membrane receptors where you have
> to deal with separate volumes.
> 
> (His solution to the problem was change the rate constant.
> I'll ask if I can forward his reply around.)
> 
> What may be a bit confusing is the rate law I used,
> mass-action, has assumptions about molecular interactions.
> If you remove that, perhaps it won't require thinking that
> the reactants have to be in the same space.
> 
> More generally, I think you would agree that modelers should
> be able to express a compound rate expression that involves
> whatever compartmental organization they feel is appropriate
> for their problem.  (And indeed, Nicolas will attest to the
> fact that models in the wild do this.)
> 
> MH
> 

Re: [SOSlib-devel] Request for help with a multicompartment problem

[Rainer M. <ra...@tb...>]

Hi Michael,

Michael Hucka wrote:


>... and goes on to give an example of extracellular
>signal molecule binding to membrane receptors where you have
>to deal with separate volumes.
>
>(His solution to the problem was change the rate constant.
>I'll ask if I can forward his reply around.)
>

I haven't read through all the messages yet, but we ran across this
problem  when trying to translate the Schoeberl et al. 2001 model to
SBML. Tell me if this initial comment already helps you, otherwise I
will try to respond in more detail later.

d[A]/dt = -k * [A] * [B]

If A and B are ligand and receptor, and thus A is soluble and 
extracellular
while B is in a membrane compartment, the reaction itself still happens in
A's compartment!! B's reacting domain is also in this compartment! If
in the rest of the model B is described in concentrations in membrane or
cytosol, the modeller will have to convert B's concentration to 
extracellular by a volume ratio Vb/Va and the equation would read


d[A]/dt = -k * [A] * Vb/Va * [B]

... and I think it's the modeller's job to write the kinetic law this way.

I agree with Andrew that in general a reaction between species in 
different compartments makes no physical sense and therefor there can also 
not be correct rate law for such a situation.

Sorry, if this misses the topic our your question and is not really 
helpful or if the above example is also wrong and already covered in your 
previous messages. I'll go through the messages later. But let me know if 
it helps.

Rainer

Re: [SOSlib-devel] Request for help with a multicompartment problem

[Michael H. <mh...@ca...>]

Can I send around your replies to other people?  I mailed
some other people about the same question and would like to
construct a compilation of everyone's replies then mail them
to everyone again.

The replies are interesting in their divergence...

MH