Global Identifiability of Differential Models

In this tutorial, let us cover an example problem of querying the ODE for globally identifiable parameters.

Input System

Let us consider the following four-dimensional model with two outputs:

\[\begin{cases}x'(t) = lm - d \, x(t) - \beta \, x(t) \, v(t),\\ y'(t) = \beta \, x(t) \, v(t) - a \, y(t),\\ v'(t) = k \, y(t) - u \, v(t),\\ w'(t) = c \, x(t) \, y(t) \, w(t) - c \, q \, y(t) \, w(t) - b \, w(t),\\ z'(t) = c \, q \, y(t) \, w(t) - h \, z(t),\\ y_1(t) = w(t),\\ y_2(t) = z(t)\end{cases}\]

This model describes HIV dynamics[1]. Let us run a global identifiability check on this model to get the result with probability of correctness being p=0.99. To do this, we will use assess_identifiability(ode, p) function.

Global identifiability needs information about local identifiability first, hence the function we chose here will take care of that extra step for us.

julia> using StructuralIdentifiability
julia> ode = @ODEmodel( x'(t) = lm - d * x(t) - beta * x(t) * v(t), y'(t) = beta * x(t) * v(t) - a * y(t), v'(t) = k * y(t) - u * v(t), w'(t) = c * x(t) * y(t) * w(t) - c * q * y(t) * w(t) - b * w(t), z'(t) = c * q * y(t) * w(t) - h * z(t), y1(t) = w(t), y2(t) = z(t) )[ Info: Summary of the model: [ Info: State variables: w, y, v, z, x [ Info: Parameters: b, c, h, lm, d, k, u, q, a, beta [ Info: Inputs: [ Info: Outputs: y1, y2 w'(t) = -b*w(t) + c*w(t)*x(t)*y(t) - c*w(t)*q*y(t) v'(t) = k*y(t) - v(t)*u x'(t) = lm - x(t)*d - x(t)*v(t)*beta z'(t) = c*w(t)*q*y(t) - h*z(t) y'(t) = x(t)*v(t)*beta - a*y(t) y2(t) = z(t) y1(t) = w(t)
julia> global_id = assess_identifiability(ode, 0.99)[ Info: Assessing local identifiability [ Info: Local identifiability assessed in 9.00206953 seconds [ Info: Assessing global identifiability [ Info: Computing IO-equations [ Info: Computing in lex order, result is in lex order [ Info: Using fglm: false [ Info: Selected tablesize 512 [ Info: 1: selected lucky prime 2147483647 [ Info: CRT modulo (1, 2147483647) [ Info: Reconstructing modulo 2147483647 [ Info: Heuristic check passed! [ Info: Randomized check failed. [ Info: 2: selected lucky prime 2147483659 [ Info: CRT modulo (2147483647, 2147483659) [ Info: Reconstructing modulo 4611686039902224373 [ Info: Heuristic check failed. [ Info: 3: selected lucky prime 2147483693 [ Info: CRT modulo (4611686039902224373, 2147483693) [ Info: Reconstructing modulo 9903520567925774155444649489 [ Info: 4: selected lucky prime 2147483713 [ Info: CRT modulo (9903520567925774155444649489, 2147483713) [ Info: Reconstructing modulo 21267649120981110191733715050621272657 [ Info: Heuristic check failed. [ Info: 5: selected lucky prime 2147483743 [ Info: CRT modulo (21267649120981110191733715050621272657, 2147483743) [ Info: Reconstructing modulo 45671930739135174346839766056203605080877915151 [ Info: Heuristic check failed. [ Info: 6: selected lucky prime 2147483777 [ Info: CRT modulo (45671930739135174346839766056203605080877915151, 2147483777) [ Info: Reconstructing modulo 98079730326560405919904968824172512120100095704355005327 [ Info: Heuristic check failed. [ Info: 7: selected lucky prime 2147483783 [ Info: CRT modulo (98079730326560405919904968824172512120100095704355005327, 2147483783) [ Info: 8: selected lucky prime 2147483813 [ Info: CRT modulo (210624630317301765882893117451031048172285903861850336414611112041, 2147483813) [ Info: Reconstructing modulo 452312984225514596069828623335196996110407213751397785678981819797976892333 [ Info: Heuristic check passed! [ Info: Randomized check passed! [ Info: Success! [ Info: Computing in lex order, result is in lex order [ Info: Using fglm: false [ Info: Computing in lex order, result is in lex order [ Info: Using fglm: false [ Info: Computing in lex order, result is in lex order [ Info: Using fglm: false [ Info: Computing in lex order, result is in lex order [ Info: Using fglm: false [ Info: Computing in lex order, result is in lex order [ Info: Using fglm: false ┌ Info: Computed in 16.72824108 seconds │ :ioeq_time = :ioeq_time └ ioeq_time = 16.72824108 [ Info: Computing Wronskians ┌ Info: Computed in 6.358604965 seconds │ :wrnsk_time = :wrnsk_time └ wrnsk_time = 6.358604965 [ Info: Dimensions of the wronskians [40, 18] ┌ Info: Ranks of the Wronskians computed in 0.03532328 seconds │ :rank_time = :rank_time └ rank_times = 0.03532328 [ Info: Assessing global identifiability using the coefficients of the io-equations [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Selected tablesize 65536 [ Info: 1: selected lucky prime 2147483647 [ Info: CRT modulo (1, 2147483647) [ Info: Reconstructing modulo 2147483647 [ Info: Heuristic check failed. [ Info: 2: selected lucky prime 2147483659 [ Info: CRT modulo (2147483647, 2147483659) [ Info: Reconstructing modulo 4611686039902224373 [ Info: Heuristic check failed. [ Info: 3: selected lucky prime 2147483693 [ Info: CRT modulo (4611686039902224373, 2147483693) [ Info: Reconstructing modulo 9903520567925774155444649489 [ Info: 4: selected lucky prime 2147483713 [ Info: CRT modulo (9903520567925774155444649489, 2147483713) [ Info: Reconstructing modulo 21267649120981110191733715050621272657 [ Info: Heuristic check failed. [ Info: 5: selected lucky prime 2147483743 [ Info: CRT modulo (21267649120981110191733715050621272657, 2147483743) [ Info: Reconstructing modulo 45671930739135174346839766056203605080877915151 [ Info: 6: selected lucky prime 2147483777 [ Info: CRT modulo (45671930739135174346839766056203605080877915151, 2147483777) [ Info: Reconstructing modulo 98079730326560405919904968824172512120100095704355005327 [ Info: 7: selected lucky prime 2147483783 [ Info: CRT modulo (98079730326560405919904968824172512120100095704355005327, 2147483783) [ Info: 8: selected lucky prime 2147483813 [ Info: CRT modulo (210624630317301765882893117451031048172285903861850336414611112041, 2147483813) [ Info: Reconstructing modulo 452312984225514596069828623335196996110407213751397785678981819797976892333 [ Info: Heuristic check failed. [ Info: 9: selected lucky prime 2147483857 [ Info: CRT modulo (452312984225514596069828623335196996110407213751397785678981819797976892333, 2147483857) [ Info: 10: selected lucky prime 2147483867 [ Info: CRT modulo (971334831935788242577832613368869049061991281227475155931159242212638377544144568381, 2147483867) [ Info: 11: selected lucky prime 2147483869 [ Info: CRT modulo (2085925881037261630864178029036094802896257759330662854505473835259586299281105540913875809327, 2147483869) [ Info: 12: selected lucky prime 2147483887 [ Info: CRT modulo (4479492181457132340213454847299227207974448018628682717127999033421525005319580445599067818799052246163, 2147483887) [ Info: Reconstructing modulo 9619637281621671881834996425177135896677125027724171971067756840824299427891388292524278203091200389506200075581 [ Info: Heuristic check failed. [ Info: 13: selected lucky prime 2147483893 [ Info: CRT modulo (9619637281621671881834996425177135896677125027724171971067756840824299427891388292524278203091200389506200075581, 2147483893) [ Info: 14: selected lucky prime 2147483929 [ Info: CRT modulo (20658016118784845285971654106780479009986238218584837754630069827310747864405871311604659752589335646499890885945580116833, 2147483929) [ Info: 15: selected lucky prime 2147483951 [ Info: CRT modulo (44362757620113410260375536343857928604867277085576528201140520294297636327782679774913158020198714437645340777821705269480809876857, 2147483951) [ Info: 16: selected lucky prime 2147483993 [ Info: CRT modulo (95268310011296503334035195531452119103056298026345647894248167227793970831147880232398299267103673085675359550297968805662169313032693822007, 2147483993) [ Info: 17: selected lucky prime 2147483999 [ Info: CRT modulo (204587170789420890086711714502418553819742917389414771138112095891274737061798978614956467676528769422992752228784566390572836365393516268430023633951, 2147483999) [ Info: 18: selected lucky prime 2147484007 [ Info: CRT modulo (439347675670961559937551129419800091128618245387347072993342744994866141553145580730162196416906240199037398104286443541908350538727893524799664605101605650049, 2147484007) [ Info: 19: selected lucky prime 2147484037 [ Info: CRT modulo (943492107016012944277662969173807884835850261977729359411465162345954336091178275160750699321298855245933309223847255653156617041667985469306137618419668743501066266343, 2147484037) [ Info: 20: selected lucky prime 2147484041 [ Info: CRT modulo (2026134238852383501221651721966775511189722802919441848842357150919550408286738322428905735729076025747015490724796841241411843757904162649282884001681435793496387299450782866691, 2147484041) [ Info: Reconstructing modulo 4351090942859175723685201076583819542509546642483289579016496306421962976690804699359187424571554264967420869711284739533142562778484656896803273488065100532499663316725644271175152978331 [ Info: Heuristic check passed! [ Info: Randomized check passed! [ Info: Success! [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false [ Info: Computing in degrevlex order, result is in degrevlex order [ Info: Using fglm: false ┌ Info: Computed in 1.234168944 seconds │ :check_time = :check_time └ check_time = 1.234168944 [ Info: Global identifiability assessed in 25.038252194 seconds Dict{Any, Symbol} with 10 entries: d => :globally lm => :nonidentifiable a => :globally c => :nonidentifiable k => :nonidentifiable q => :nonidentifiable h => :globally beta => :nonidentifiable b => :globally u => :globally

We also note that it's usually inexpensive to obtain the result with higher probability of correctness. For example, taking p=0.9999 in the system above will result only in a slight slowdown.

Note on the probability of correctness

Currently, the probability of correctness does not include the probability of correctness of the modular reconstruction for Groebner bases. This probability is ensured by additional check modulo a large prime and can be neglected for practical purposes. However, in the future versions, we plan to eliminate this possible error.