Identifiability with known initial conditions.

[pogudingleb]

It seems that there is an issue with the structural identifiability assessment if the initial conditions are given as specific numbers. I did the following

1. enabled using known initial conditions by setting opts.replaceICs option to be equal to one in the options.m
2. tried the following system (full matlab code below)
   S' = 1 + a * S^2
   S(0) = 0
   y = S
3. I get the result that the parameter a is not identifiable.

However, since S(0) = 0 and S'(0) = 1, one can see that
y''' = (S')'' = 2 * a * S * S'' + 2 * a * S' * S'
being evaluated at t = 0 yields y'''(0) = 2a, so a is globally identifiable.

It seems that the issue is that, in the case of known initial conditions, one should compute the observability matrix further than in the generic case.

==========================

clear all;

% one state
syms S
x = [S].';

% 1 parameter:
syms a
p = [a].';

% 1 output:
h = [S].';

% no inputs
u = [];
w = [];

% dynamic equations:
f = [1 + a * S^2];

% initial conditions:
ics  = [0.];   

% which initial conditions are known:
known_ics = [1]; 

save('example','x','p','u','w','h','f','ics','known_ics');

[afvillaverde]

Thanks for warning us of this issue! Indeed, this example shows that, for specific values of the initial conditions, the OIC may wrongly classify a parameter as unidentifiable. The cause of the error is that for some specific initial conditions it may be necessary to build the observability-identifiability matrix with more Lie derivatives than for the generic case. We were not aware of this problem until now. It seems to be a rare occurrence, but of course we must take it into account.

To solve this issue we need to determine the upper bound on the number of Lie derivatives that have to be calculated. For the generic case, there is a theoretical results that establishes that the bound equals the dimension of the model, i.e. the total number of unknowns. For particular cases we do not have an answer yet.

In summary: until we fix this issue, the results obtained with opts.replaceICs = 1 should only be trusted if they report that the model is FISPO. A result reporting that the model is FISPO with with opts.replaceICs = 0, and it is not FISPO with opts.replaceICs = 1, should be considered non-conclusive instead. If we do not find a solution in the short term, we will modify the output of the toolbox to reflect this fact.

[afvillaverde]

We have just updated the toolbox, removing a number of options. In the new version it is no longer possible to replace ICs with specific numerical conditions, so this issue does not apply any more.