Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Cannot increase weighting of likelihood in examples/Shin_PLOS2019 #11

Closed
paulflang opened this issue Oct 4, 2020 · 6 comments
Closed
Assignees

Comments

@paulflang
Copy link
Owner

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
Copy link
Owner Author

This issue refers to commit ff8dfac

@sshin23
Copy link
Collaborator

sshin23 commented Oct 25, 2020

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
Copy link
Owner Author

paulflang commented Oct 25, 2020

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:
image

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
Copy link
Owner Author

I am trying MA97 now.

@paulflang
Copy link
Owner Author

MA97 exits with
EXIT: Maximum Number of Iterations Exceeded.

@paulflang
Copy link
Owner Author

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 07f2a73 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).

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants