@@ -947,7 +947,133 @@ ode = @ODEmodel(x'(t) = alpha * x(t)^2, y(t) = x(t)^2)
947947ident_funcs = [alpha^ 2 , x * alpha]
948948push! (test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = true ))
949949
950- # TODO : verify that Maple returns the same
950+ # Examples from the paper "An algorithm for finding globally identifiable parameter combinations of nonlinear ODE models using Gröbner Bases"
951+ # http://dx.doi.org/10.1016/j.mbs.2009.08.010
952+
953+ # Section 6
954+ ode = @ODEmodel (
955+ x1' (t) = - (k01 + k21) * x1 (t) + k12 * x2 (t) + u (t),
956+ x2' (t) = k21 * x1 (t) - (k02 + k12) * x2 (t),
957+ y (t) = x1 (t) / v
958+ )
959+ ident_funcs = [v, k02 + k12, k01 + k21, k12 * k21]
960+ ident_funcs_states = [v, x1, k02 + k12, k01 + k21, k12 * k21, x2 * k12]
961+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = false ))
962+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs_states, with_states = true ))
963+
964+ # Example 1
965+ ode = @ODEmodel (
966+ x1' (t) = - (k21 + V_M / (K_M + x1 (t))) * x1 (t) + k12 * x2 (t) + b1 * u (t),
967+ x2' (t) = k21 * x1 (t) - (k02 + k12) * x2 (t),
968+ y (t) = c * x1 (t)
969+ )
970+ ident_funcs = [k21, k12, k02, b1 * c, V_M * c, K_M * c]
971+ ident_funcs_states = [k21, k12, k02, b1 * c, V_M * c, K_M * c, x2 * c, x1 * c]
972+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = false ))
973+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs_states, with_states = true ))
974+
975+ # Example 2
976+ ode = @ODEmodel (
977+ x1' (t) = k13 * x3 (t) + k12 * x2 (t) - (k12 + k31) * x1 (t) + u (t),
978+ x2' (t) = k21 * x1 (t) - (k12 + k02) * x2 (t),
979+ x3' (t) = k31 * x1 (t) - (k13 + k03) * x3 (t),
980+ y (t) = x1 (t) / v
981+ )
982+ ident_funcs = [
983+ v,
984+ k12 + k31,
985+ k02 + k03 + k13 - k31,
986+ k12* k21 + k13* k31,
987+ k02* k03 + k02* k13 + k03* k12 + k12* k13,
988+ k02* k13* k31 + k03* k12* k21 + k12* k13* k21 + k12* k13* k31,
989+ ]
990+ ident_funcs_states = [
991+ v,
992+ x1,
993+ k12 + k31,
994+ k02 + k03 + k13 - k31,
995+ k12* k21 + k13* k31,
996+ x2* k12 + x3* k13,
997+ k02* k03 + k02* k13 + k03* k12 + k12* k13,
998+ k02* k13* k31 + k03* k12* k21 + k12* k13* k21 + k12* k13* k31,
999+ x2* k03* k12 + x2* k12* k13 + x3* k02* k13 + x3* k12* k13,
1000+ ]
1001+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = false ))
1002+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs_states, with_states = true ))
1003+
1004+ # Example 3
1005+ ode = @ODEmodel (
1006+ x1' (t) = - a31 * x1 (t) + a13 * x3 (t) + u (t),
1007+ x2' (t) = - a42 * x2 (t) + a24 * x4 (t),
1008+ x3' (t) = a31 * x1 (t) - (a03 + a13 + a43) * x3 (t),
1009+ x4' (t) = a42 * x2 (t) + a43 * x3 (t) - (a04 + a24) * x4 (t),
1010+ y1 (t) = x1 (t),
1011+ y2 (t) = x2 (t)
1012+ )
1013+ ident_funcs = [a31, a13, a03 + a43, a04 + a24 + a42, a24* a43, a04* a42]
1014+ ident_funcs_states =
1015+ [a31, a13, x3, x2, x1, a03 + a43, a04 + a24 + a42, a24* a43, a04* a42, x2* a42 - x4* a24]
1016+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = false ))
1017+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs_states, with_states = true ))
1018+
1019+ # Example 4
1020+ ode = @ODEmodel (
1021+ S' (t) = mu * N - S * (mu + beta / N * I (t)),
1022+ I' (t) = I (t) * (beta / N * S (t) - mu - nu),
1023+ y (t) = k * I (t)
1024+ )
1025+ ident_funcs = [nu, mu, beta, N* k]
1026+ ident_funcs_states = [nu, mu, beta, N* k, I* k, S* k]
1027+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = false ))
1028+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs_states, with_states = true ))
1029+
1030+ # Example 5
1031+ ode = @ODEmodel (
1032+ x1' (t) = p1 * x1 (t) - p2 * x1 (t) * x2 (t),
1033+ x2' (t) = p3 * x2 (t) * (1 - p4 * x2 (t)) + p5 * x1 (t) * x2 (t),
1034+ y (t) = x1 (t)
1035+ )
1036+ ident_funcs = [p5, p3, p1, p2// p4]
1037+ ident_funcs_states = [p5, p3, p1, x1, x2* p4, p2// p4]
1038+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = false ))
1039+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs_states, with_states = true ))
1040+
1041+ # Example 6
1042+ ode = @ODEmodel (
1043+ x1' (t) = - beta * x1 (t) * x4 (t) - d * x1 (t) + s,
1044+ x2' (t) = beta * q1 * x1 (t) * x4 (t) - k1 * x2 (t) - mu1 * x2 (t),
1045+ x3' (t) = beta * q2 * x1 (t) * x4 (t) + k1 * x2 (t) - mu2 * x3 (t),
1046+ x4' (t) = - c * x4 (t) + k2 * x3 (t),
1047+ y1 (t) = x1 (t),
1048+ y2 (t) = x4 (t)
1049+ )
1050+ ident_funcs = [
1051+ s,
1052+ d,
1053+ beta,
1054+ c + k1 + mu1 + mu2,
1055+ k2* q2,
1056+ c* k1 + c* mu1 + c* mu2 + k1* mu2 + mu1* mu2,
1057+ c* k1* mu2 + c* mu1* mu2,
1058+ k1* k2* q1 + k1* k2* q2 + k2* mu1* q2,
1059+ ]
1060+ ident_funcs_states = [
1061+ s,
1062+ d,
1063+ beta,
1064+ x4,
1065+ x1,
1066+ c + k1 + mu1 + mu2,
1067+ k2* q2,
1068+ x3* k2 - x4* c,
1069+ c* k1 + c* mu1 + c* mu2 + k1* mu2 + mu1* mu2,
1070+ c* k1* mu2 + c* mu1* mu2,
1071+ k1* k2* q1 + k1* k2* q2 + k2* mu1* q2,
1072+ x2* k1* k2 + x3* k1* k2 + x3* k2* mu1 + x4* k1* mu2 + x4* mu1* mu2,
1073+ ]
1074+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = false ))
1075+ push! (test_cases, (ode = ode, ident_funcs = ident_funcs_states, with_states = true ))
1076+
9511077@testset " Identifiable functions of parameters" begin
9521078 p = 0.99
9531079 for case in test_cases
0 commit comments