*very* long identifiability runtime for relatively small model

[TorkelE]

I have the following model:

ode = @ODEmodel(
    X'(t) = p - X(t),
    Y'(t) = -Y(t) + ((K1^3)*v1) / (K1^3 + X(t)^3),
    Z'(t) = -Z(t) + (v2*(Y(t)^3)*(X(t)^3)) / ((K2^3 + X(t)^3)*(K3^3 + Y(t)^3)),
    y(t) = Z(t)
)

for which the runtime of assess_identifiability is incredibly long. I.e.

assess_identifiability(ode)

takes at least 24 hours (I have not actually managed to complete it).

The model is essentially a incoherent feedforward loop (X deactives Y and activets Z, Y activates Z). In Catalyst it can be implemented like

rn = @reaction_network begin
    (p,1.0), 0 <--> X
    hillr(X, v1, K1, 3), 0 --> Y
    v2*hill(X, 1.0, K2, 3)*hill(Y, 1.0, K3, 3), 0 --> Z
    1.0, (Y,Z) --> 0
end

[pogudingleb]

Hi @TorkelE , thanks a lot for reporting!
I have checked briefly, this seemingly innocent model indeed triggers complicated computations (because of high nonlinearity appearing in denominators), so this is in some sense "correct behavior".
But, of course, this is annoying (if not embarrasing) to take so long on such a small model, I will think about this, there should be a way to get this work.

[TorkelE]

Thanks for the reply! Just happy to know that the sluggishness is a property of the model and nothing else.

Generally, if I have something which takes like a day, is there hope in putting it on an HPC and let run for a week, or at that point it is probably better to give up?

[pogudingleb]

Yes, in some cases passing to a machine with more memory and waiting longer may allow to get ultimately result. Looking at the internal debugging info in this example, I am skeptical about this particular case.

However! I have found a workaround (but not sure how to automatize it). From the first equation, we know that p - X(t) is proportional to $e^{-t}$. So we can take any known output $y_2(t) = b e^{-t}$. Then, for some (unknown) constant a, we will have y2(t) = a(p - X(t)). This will give us the following system (equivalent identifiability-wise to the original one):

ode = @ODEmodel(
    X'(t) = p - X(t),
    Y'(t) = -Y(t) + ((K1^3)*v1) / (K1^3 + X(t)^3),
    Z'(t) = -Z(t) + (v2*(Y(t)^3)*(X(t)^3)) / ((K2^3 + X(t)^3)*(K3^3 + Y(t)^3)),
    y(t) = Z(t),
    y2(t) = a * (p - X(t)),
)

After this rewriting, the identifiability assessment is very quick.

PS. This seem to be also related to first integrals as we essentially use that $(p - X(t)) e^t$ is a first integral of the model.

[TorkelE]

Thanks a lot for that extra!

Works flawlessly, could compute identifiability for my model, and a couple of variants of it, without any problems.

[pogudingleb]

A bit longer explanation/justification for myself, in order not to forget.

1. Augment the system with an observed state satisfying W'(t) = -W(t):

ode = @ODEmodel(
    X'(t) = p - X(t),
    Y'(t) = -Y(t) + ((K1^3)*v1) / (K1^3 + X(t)^3),
    Z'(t) = -Z(t) + (v2*(Y(t)^3)*(X(t)^3)) / ((K2^3 + X(t)^3)*(K3^3 + Y(t)^3)),
    y(t) = Z(t),
    W'(t) = -W(t),
    y2(t) = W(t),
)

The subsystem consisting of the last two equations is not coupled with the original, so this transformation does not change identifiability.

2. Remark that the new system has a first integral C = W(t) / (p - X(t)). Use it to eliminate W from the original model:

ode = @ODEmodel(
    X'(t) = p - X(t),
    Y'(t) = -Y(t) + ((K1^3)*v1) / (K1^3 + X(t)^3),
    Z'(t) = -Z(t) + (v2*(Y(t)^3)*(X(t)^3)) / ((K2^3 + X(t)^3)*(K3^3 + Y(t)^3)),
    y(t) = Z(t),
    y2(t) = C * (p - X(t)),
)

So we get the extended model as described in the post above.