Cannot increase weighting of likelihood in `examples/Shin_PLOS2019`

[paulflang]

In the petab_suite branch, when I reduce the measurement error​ (examples/Shin_PLOS2019/julia_code.jl line 252) from sigma = 1.0 to 0.1, 0.5 or 0.8, JuMP either exits with Restoration failed! or Maximum number of iterations exceeded. Any ideas how to solve that?

[paulflang]

This issue refers to commit ff8dfac9cc935a9abd98fae949b2b5b7bb26d4ea

[sshin23]

Seems to work if the linear solver is changed to ma27 or ma97. This is a bit mysterious behavior, but sometimes using a different linear solver fixes the issue.

[paulflang]

Thanks for the hint. I tried MA27, and it works without restoration error, but the estimates are not great. See for instance this image for the growth of species one alone:

Here is the code I used to run it:

import os
import pandas as pd
import pickle
import pkg_resources
import re
import shutil
import tempfile
import unittest
from SBML2JuliaMP import core
import importlib
importlib.reload(core)

FIXTURES = os.path.join('/media/sf_DPhil_Project/Project07_Parameter Fitting/df_software',
    'SBML2JuliaMP', 'examples')
PETAB_YAML = os.path.join(FIXTURES, 'Shin_PLOS2019', 'Shin_PLOS2019.yaml')

problem = core.SBML2JuliaMPProblem(PETAB_YAML, optimizer_options={'linear_solver': 'MA27'}, n_starts=1, t_steps=48)
problem.insert_custom_code({'    # Model initial assignments': '    # Custom code\n    bindings = [Int(i) for i in range(13, stop=210, length=198) if mod(i,3) != 0]\n    @constraint(m, [j in bindings], A[j, 3] == 20*A[j+1, 1])\n    @constraint(m, [j in bindings], B[j, 3] == 20*B[j+1, 1])\n\n    # Model initial assignments'})
problem.write_jl_file()
problem.optimize()
problem.plot_results('c_1', path='plot.pdf')
problem.write_results()
problem.write_optimized_parameter_table()

Here is the generated jl_file: julia_code.txt

When fitting, I noticed the following messages:

 818  1.5866191e+03 1.07e+00 2.22e+01  -1.7 2.83e+01    -  2.87e-01 1.00e+00f  1
MA27BD returned iflag=-4 and requires more memory.
 Increase liw from 7488620 to 14977240 and la from 8899304 to 31029470 and factorize again.
MA27BD returned iflag=-4 and requires more memory.
 Increase liw from 14977240 to 29954480 and la from 31029470 to 125264652 and factorize again.
WARNING: Problem in step computation; switching to emergency mode.
 819r 1.5866191e+03 1.07e+00 1.00e+03   0.0 0.00e+00  19.4 0.00e+00 0.00e+00R  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
 820r 1.5885163e+03 1.32e-02 9.63e+02   0.0 1.09e+03    -  7.31e-02 1.02e-03f  1
 821  1.5919101e+03 2.69e-01 3.62e+03  -1.7 5.46e+00    -  8.00e-01 9.99e-01h  1
 822  1.5959770e+03 4.63e-02 1.31e+00  -1.7 3.43e+00    -  1.00e+00 1.00e+00f  1
 823  1.5432241e+03 2.77e+00 2.80e+04  -2.5 4.42e+01    -  6.66e-01 6.84e-01f  1
 824  1.5068525e+03 2.83e+00 1.08e+05  -2.5 2.85e+01    -  4.80e-01 6.71e-01f  1
 825  1.4822980e+03 2.38e+00 1.03e+05  -2.5 1.97e+01    -  3.88e-01 6.22e-01f  1
 826  1.4731651e+03 1.27e-01 1.58e+05  -2.5 7.88e-01    -  3.41e-02 1.00e+00f  1
 827  1.4818965e+03 6.17e-01 2.64e+04  -2.5 6.34e+00    -  8.34e-01 1.00e+00h  1
Scaling factors are invalid - setting them all to 1.
 828  1.4768289e+03 2.21e-01 1.02e+00  -2.5 2.92e+00    -  1.00e+00 1.00e+00h  1
Scaling factors are invalid - setting them all to 1.

[paulflang]

I am trying MA97 now.

[paulflang]

MA97 exits with EXIT: Maximum Number of Iterations Exceeded.

[paulflang]

Works now. Previously it appeared that either the weights on the priors compared to the likelihood was too strong, or the restoration failed. I do not know what exactly went wrong (maybe an unfortunate combination of solvers and priors), but in commit 07f2a73eceaf41756aa9be27ce89e83c0e7cee6e everything works fine now. The fits actually look pretty good, but this because I set the measurement noise quite low (0.05) and the standard deviation of the parameter prior quite high (3).